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SECTION  I 


1.1  Objectives 

Tnis  report  presents  a  series  of  investigations  that  have  been  carried  out  for  the 
application  of  a  packed  bed  (with  encapsulated  phase  change  material-PCM)  as  an  energy 
storage/release  unit  for  spacecraft  thermal  management  systems.  In  this  regard,  relevant 
studies  performed  in  the  Department  of  Mechanical  Engineering  at  The  Ohio  State 
University  before  and  during  the  course  of  the  Contract  F33615-89-C-2949  have  been 
assembled  in  this  report  for  integrity  of  the  material  presented.  The  fundamental  studies 
performed  form  an  integral  part  of  the  preliminary  investigations  for  the  actual  application, 
and  these  will  be  presented  sequentially  in  the  present  report. 

A  set  of  investigations  related  to  the  single  and  multiphase  transport  processes  in 
packed  beds,  and  the  behavior  of  sensible  and  latent  heat  storage  packed  beds  was  carried 
out.  The  main  objective  of  this  work  was  to  analyze  different  aspects  of  the  applications  of 
packed  beds  by  considering  a  set  of  problems  which  were  designed  to  provide  physical 
insight  for  a  number  of  outstanding  points  related  to  transport  phenomena  in  porous  media 
as  well  as  different  packed  bed  energy  storage  systems  including  sensible  heat  storage  and 
latent  heat  storage  ones. 

1.2  ln(r9.d»Cti.9n 

Packed  beds  with  single-phase  flow  have  been  analyzed  in  applications  ranging  from 
chemical  catalytic  reactors  and  pebble  bed  or  rock  pile  heat  storage  units  to  fixed-bed 
nuclear  propulsion  systems  and  spacecraft  thermal  management  systems  employing  packed 
beds  of  encapsulated  phase  change  materials.  A  major  part  of  the  studies  conducted  to  date 
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on  packed  bed  applications,  w  ith  single-phase  flow  concentrates  on  utilizing  incompressible 
fluids,  liquid  or  gas,  as  the  heat  transporting  medium.  In  these  investigations  a  constant 
mass  flow  rate  is  assumed  at  every  cross  section  of  the  packed  bed  and,  therefore,  there  is 
no  need  for  solving  any  continuity  or  momentum  equations.  The  problem  reduces  to 
solving  the  governing  energy  equations  for  the  solid  and  fluid  phases.  For  this  category  of 
problems,  the  two  commonly  employed  models  are  the  so  called  single-phase  conductivity 
and  two-phase  Schumann  models.  The  main  distinction  between  these  two  models  is  that, 
where  local  thermal  equilibrium  is  assumed  to  prevail  between  the  solid  and  the  fluid 
phases  in  the  single-phase  model,  no  such  assumption  is  made  in  the  two-phase  model. 
Therefore,  the  single  phase  conductivity  model  reduces  to  one  governing  energy  equation 
in  which  conduction  in  both  phases  is  taken  care  of  by  the  use  of  an  effective  thermal 
conductivity,  where  in  the  two-phase  model  there  are  two  governing  energy  equations, 
each  of  them  possessing  a  fluid-to-solid  convective  heat  transfer  term.  The  origin  of  the 
two-phase  models  is  the  classical  Schumann  model  (Schumann,  1929).  Riaz  (1977) 
presented  an  investigation  and  comparison  of  both  the  single-phase  conductivity  model  as 
well  as  the  two-phase  Schumann  model.  From  a  theoretical  point  of  view,  both  of  these 
models  are  quite  inadequate  because  they  incorporate  many  simplifying  assumptions,  such 
as  neglecting  the  transient  heat  storage  and  the  conduction  terms  in  the  fluid  phase. 

The  second  major  class  of  problems  dealing  with  the  single-phase  flow  through 
porous  media  studied  to  date  concentrates  on  the  compressible  gas  flow  through  a  porous 
medium.  It  is  essential  to  consider  a  compressible  working  fluid-hence,  an  equation  of 
state  for  the  vapor  phase-in  studying  multiphase  flow  with  phase  change  in  porous  media. 
In  this  category  of  problcms  too,  there  is  a  lack  of  analysis  of  the  complete  transient,  non- 
isothermal  and  nonthermal  equilibrium  flow  of  a  compressible  gas.  One-dimensional, 
transient,  isothermal  flow  of  an  ideal  gas  through  a  porous  medium  has  been  studied 
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analytically  by  the  use  of  perturbation  methods  by  Kidder  (1957).  Morrison  (1972)  also 
obtained  analytical  solutions  by  perturbation  methods  for  one-dimensional  isothermal  and 
adiabatic  flows  of  an  ideal  gas.  In  both  of  these  investigations  Darcy's  formulation  has 
been  used  for  fluid  flow.  Nilson  (1981),  on  the  other  hand,  obtained  an  analytical  solution 
for  one-dimensional  isothermal  flow  through  a  porous  medium  in  a  study  which  utilized 
Ergun’s  equation  (Ergun,  1952)  to  account  for  the  inertial  effects.  In  all  of  these 
investigations  it  was  assumed  that  the  gas  phase  was  in  local  thermal  equilibrium  with  the 
solid  phase.  Goldstein  and  Siegel  (1971)  performed  :n  analytical  investigation  on  a  steady, 
nonisothermal  ideal  gas  flow  in  a  porous  medium.  They  also  utilized  the  local  thermal 
equilibrium  assumption  in  the  governing  energy  equation  and  Darcy's  formulation  for  fluid 
flow.  Most  of  these  studies  which  incorporated  Darcy's  formulation  for  fluid  motion  have 
the  drawback  of  not  being  applicable  to  high  speed  flows  in  which  the  inertia  effects  have 
to  be  accounted  for. 

It  can  be  seen  that  the  models  used  in  both  categories  of  problems  described  above  are 
incomplete  in  one  aspect  or  another  since  they  concentrate  on  either  the  momentum 
transport  in  an  isothermal  or  adiabatic  flow  or  the  energy  transport  in  an  incompressible 
flow.  In  addition  they  incorporate  some  simplifying  assumptions  such  as  local  thermal 
equilibrium  between  the  solid  and  the  fluid  phases,  neglect  of  conduction  and  transient 
storage  terms  in  the  gas  phase  energy  equation,  and  neglect  of  inertial  effects  in  the  fluid 
phase  momentum  equations.  Therefore,  their  application  is  limited  to  certain  special  cases. 
These  cause  a  need  for  a  more  rigorous  formulation  of  the  transport  processes  in  a  porous 
medium  for  establishing  a  more  flexible  and  reliable  model  with  a  wide  range  of 
applicability.  The  main  objective  of  Section  II  is  to  present  a  rigorous  analysis  and 
investigation  for  transient  transport  processes  in  a  compressible  gas  flow  through  a  packed 
bed  with  no  assumption  of  local  thermal  equilibrium  between  the  solid  and  fluid  phases. 
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Ergun-Forchheimer  relation  will  be  employed  as  the  vapor  phase  momentum  equation  in 
order  to  account  for  the  inertia  effects  in  addition  to  the  viscous  effects.  This  is  very 
essential  for  non-Darcy  regime  flows  in  porous  media.  More  specifically,  it  is  essential  in 
flows  with  particle  Reynolds  numbers  greater  than  0.1  in  which  Darcy  formulation 
becomes  inaccurate.  This  point  has  been  neglected  in  many  of  the  previous  studies  dealing 
with  single  and  multiphase  transport  phenomena  in  porous  media.  Likewise,  the  majority 
of  the  previous  works  on  multiphase  transport  in  porous  media  almost  exclusively  employ 
local  thermal  equilibrium  (LTE)  between  the  solid  and  fluid  phases  considered.  This 
assumption  is  usually  not  satisfactory  for  the  step  change  problems  in  which,  during  the 
early  stages  of  the  transport  processes,  there  may  be  considerable  difference  between  the 
temperatures  of  the  flowing  fluid  and  solid  particles.  This  is  also  true  even  during  the  later 
stages  of  the  transport  processes  in  high-speed  flows  or  high  permeability  porous  media  in 
which  the  fluid-to-solid  interaction  time  or  interaction  surface  area  respectively  may  not  be 
large  enough  for  the  temperatures  of  the  fluid  and  solid  phases  to  become  close  enough  for 
LTE  to  be  a  reasonable  assumption.  The  thermal  interactions  between  the  vapor  and  the 
solid  phases,  and  how  the  fluid  flow  as  well  as  the  pressure  and  density  fields  are  affected 
by  these  interactions  will  be  analyzed  in  detail.  Another  objective  of  Section  II  is  to  gain  a 
better  understanding  of  the  situations  under  which  the  local  thermal  equilibrium  assumption 
and  one-  or  two-dimensional  consideration  of  certain  field  variables  would  be  justifiable. 
In  this  regard  the  effect  of  certain  characteristic  nondimensional  parameters,  such  as  the 
particle  Reynolds  number,  the  Darcy  number  and  certain  thermophysical  parameters,  on  the 
general  qualitative  behavior  of  the  transport  phenomena  in  packed  beds  will  be  investigated 
thoroughly. 

One  further  step  into  these  studies  leads  to  multiphase  flow  of  a  single  pure  substance 
in  liquid  and  vapor  form,  and  the  accompanying  transport  processes  in  packed  beds.  The 
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flow  of  a  mixture  of  a  vapor  and  noncondensible  gases,  such  as  air  and  water  vapor, 
through  a  porous  medium,  and  the  accompanying  phase  change  (condensation)  and  the 
multiphase  transport  processes  involved  have  received  considerable  attention  in  a  number 
of  investigations  related  to  different  applications.  These  include  phase  change  in  building 
insulation  materials  (Vafai  and  Whitaker,  1986;  Vafai  and  Sarkar,  1986;  Ogniewicz  and 
Tien,  1981),  heat  pipe  technology  (Udell,  1985),  drying  of  different  porous  materials 
(Berger  and  Pei,  1973;  Plumb  et  al.,  1985),  and  phase  change  in  porous  media  (Eckert  and 
Pfender,  1980;  Motakef  and  El-Masri,  1986).  The  condensing  flow  of  a  single  vapor 
through  a  porous  medium,  on  the  other  hand,  received  relatively  little  attention  (Nilson  and 
Montoya,  1980).  A  rigorous  model  which  includes  the  basic  thermodynamics  of  the 
condensation  process  and  the  concept  of  nonthermal  equilibrium  between  the  solid  and  the 
fluid  phases  under  condensing  conditions  seem  to  be  completely  absent.  This  fact  is  the 
main  motivation  for  the  investigation  in  Section  III,  which  is  aimed  at  analyzing  the  forced 
convective  condensing  flow  of  a  vapor  through  a  packed  bed.  This  is  an  essential 
fundamental  study  for  the  application  of  thermal  energy  storage/release  system  for  rejection 
of  heat  in  pulsed  space  power  supplies  and  some  conceptual  spacecraft  thermal 
management  systems.  The  operating  conditions  of  the  packed  beds  for  such  applications 
may  require  gas/vapor  flow  at  high  speeds  as  well  as  high  pressures,  which  in  turn  will 
dictate  condensation  of  the  working  fluid.  This  may  actually  be  desirable  in  order  to 
enhance  the  amount  of  thermal  energy  stored  in  the  packed  bed. 

For  analyzing  the  problem  in  Section  III,  the  model  established  for  single-phase 
vapor  flow  in  Section  II  will  be  expanded  such  that  it  will  accommodate  the  liquid  phase  in 
addition  to  the  vapor  phase  of  the  working  fluid.  This  will  be  carried  out  by  additional 
conservation  equations  for  the  liquid  phase  as  well  as  accommodations  for  the  phase 
change  terms  in  the  already  available  conservation  equations  for  the  vapor  phase. 
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Moreover,  the  thermodynamic  aspects  of  the  condensation  process  which  will  introduce  an 
additional  coupling  relation  will  be  presented. 

A  further  stage  of  this  study  is  concerned  with  the  application  of  an  encapsulated 
phase  change  material  (PCM)  as  an  energy  storage  medium  in  packed  beds.  This  problem 
will  be  taken  up  in  Section  IV.  The  earlier  forms  of  the  packed  bed  energy  storage  units 
solely  relied  on  the  sensible  heat  capacity  of  the  solid  bed  particles  for  storing  thermal 
energy.  This  form  has  been  satisfactorily  employed  for  various  applications.  However, 
certain  applications  may  impose  a  limitation  on  the  size  and  the  weight  of  the  packed  bed 
system  utilized.  For  instance,  the  present  application  of  a  heat  rejection  system  in  pulsed 
space  power  supplies  which  incorporate  packed  beds,  the  reduction  of  the  mass  and 
volume  is  of  utmost  importance.  In  such  cases,  utilization  of  only  the  sensible  heat 
capacity  of  a  certain  material  for  energy  storage  will  most  likely  prove  to  be  inefficient.  The 
remedy  to  this  can  be  found  in  the  utilization  of  latent  heat  in  the  process  of  energy  storage. 
Recently,  encapsulated  phase  change  materials  (PCM)  have  received  considerable  attention 
as  energy  storage  materials.  The  use  of  an  encapsulated  PCM  is  very  appealing  since  it 
makes  the  utilization  of  latent  heat  storage  capacity  possible.  This  is  achieved  by  using  a 
PCM  which  has  a  melting  temperature  within  the  temperature  range  of  operation  of  the 
system  incorporating  the  packed  bed.  The  principal  advantage  of  PCMs  in  packed  beds  is 
that  the  energy  storage  density  of  the  bed  is  increased  significantly  and  thus,  the  size  and 
mass  of  the  storage  system  required  for  a  particular  application  are  reduced  proportionally. 

Different  PCMs  have  been  considered  for  use  in  the  packed  bed  energy  storage  units 
in  different  applications.  For  applications  over  450°C  significant  consideration  was  given 
to  salts  by  Marianowski  and  Maru  t.1977).  The  physical  properties  of  different  PCMs 
considered  for  storage  of  solar  cnerg\  have  been  presented  by  Lane  (1986).  A  number  of 
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studies  have  been  carried  out  on  the  analysis  of  latent  heat  storage  packed  beds  by  different 
researchers.  Anarthanarayanan  Ptal.  (1987)  investigated  the  dynamic  behavior  of  a  packed 
bed  which  utilized  encapsulated  Al-Si  shois  which  have  a  melting  temperature  of  577°C. 
Air  was  used  as  the  energy  transporting  fluid  in  their  study.  Pitts  and  Hong  (1987) 
presented  another  study  on  transient  thermal  behavior  of  a  latent  heat  storage  packed  bed 
which  utilized  an  inorganic  compound  hydrate  PCM,  namely  Na2HP03.12H20,  as  the 

heat  storage  medium.  Torab  and  Chang  (1988)  investigated  the  use  of  encapsulated  phase 
change  materials  for  thermal  energy  storage  units  in  space  power  systems.  They  reported 
an  analysis  of  a  latent  heat  storage  packed  bed  which  utilized  lithium  hydride  as  the  PCM, 
and  lithium  as  the  transport  fluid. 

In  the  majority  of  models  employed  in  the  analysis  of  this  category  of  problems,  the 
superficial  velocity  of  the  working  fluid  is  assumed  to  be  constant.  This  reduces  the 
system  of  governing  equations  to  a  set  of  energy  equations  for  the  working  fluid  and  for 
the  PCM  respectively.  As  pointed  out  before,  although  this  is  a  satisfactory  approach  when 
the  working  fluid  is  incompressible,  it  is  not  so  when  the  working  fluid  is  a  gas  or  vapor 
under  high  pressures.  Moreover,  when  the  working  fluid  itself  undergoes  phase  change 
(condensation),  a  rigorous  model  which  consists  of  the  governing  energy,  mass  and 
momentum  balance  equations  in  addition  to  the  relevant  coupling  thermodynamic  relations 
has  to  be  employed  for  analyzing  any  phase  change  and  the  corresponding  transport 
phenomena.  In  Section  IV  the  analysis  of  such  a  problem  with  a  rigorous  model  which 
will  basically  be  the  expansion  of  the  mod  J  developed  in  Section  III  will  be  presented. 
The  physical  aspects  of  the  phase  change  process  in  the  PCM  will  be  modeled  and 
analyzed. 
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The  main  objective  of  Section  IV  is  to  analyze  the  energy  storage,  characteristics  of 
packed  beds  which  utilize  phase  change  materials  in  the  bed  particles.  Time  history  of  the 
crucial  field  variables  such  as  the  temperature  profiles  of  the  working  fluid  and  the  bed 
particles,  arc  the  velocity,  density  and  pressure  of  the  working  fluid  will  be  determined. 
The  determination  of  these  quantities  is  important  since  the  amo  mt  of  energy  flowing  into 
and  out  of  the  packed  bed,  and  hence  the  amount  of  energy  stored  in  the  packed  bed  as  a 
function  of  time,  can  be  determined  from  the  t  'me  histories  of  these  variables. 

The  operation  of  the  packed  bed  as  an  energy  storage/release  system  requires  the 
analysis  of  the  transient  behavior  of  the  packed  bed  in  thermal  charging  and  thermal 
discharging  modes.  This  behavior  for  both  sensible  and  latent  heat  storage  packed  beds  is 
investigated  in  Section  V  in  order  to  determine  the  qualitative  and  quantitative  characteristics 
of  the  two  types  of  packed  beds  and  the  difference  between  their  characteristics.  This  is 
performed  by  extending  the  analysis  in  Section  IV  to  include  the  thermal  discharging  mode 
operation  of  both  types  of  packed  beds. 

The  final  section  of  the  teport  deals  with  the  analysis  of  oscillating  flow  through  a 
packed  bed  This  investigation  is  canted  out  in  order  to  determine  the  characteristics  of  the 
packed  bed  systems  under  more  realistic  boundary  conditions  than  the  ideal  conditions. 
For  this  reason,  oscillating  temperatuie  ot  piessure  inlet  boundary  conditions  are  employed 
in  the  case  stut.ies  in  Section  VI.  The  importance  of  such  boundary  conditions  on  the 
energy  storage  behavior  of  the  packed  bed  system  is  investigated. 
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SECTION  II 


ANALYSIS  OF  TRANSPORT  PHENOMENA  IN  A 
SINGLE-PHASE  FLOW  OF  A  GAS  IN  A  PACKED  BED 

In  this  section,  the  model  developed  for  analyzing  the  single-phase  fluid  flow  and 
transport  processes  in  a  packed  bed  is  presented.  The  numerical  method  of  solution  is 
discussed.  A  benchmarking  tnat  was  carried  out  for  verification  of  the  results  of  the 
numerical  code  developed  in  the  present  work  is  also  presented.  This  is  followed  by  the 
presentation  and  a  thorough  discussion  of  the  case  studies  performed  and  the  result 
obtained. 


2.1  Stfli£in£n.tJiLlh£„p.ml?km 

The  problem  that  will  be  investigated  in  this  section  is  the  transfer  and  storage  of  heat 
from  a  high  temperature  reservoir  into  a  packed  bed  which  consists  of  randomly  packed 
spherical  particles  of  uniform  size.  The  schematic  diagram  of  the  problem  considered  is 
depicted  in  Figure  2.1.  The  extent  of  the  packed  bed  in  the  z-direction  is  assumed  to  be 
large  enough  so  that  the  problem  will  essentially  be  two-dimensional.  Superheated 
Refrigerant- 12  is  used  as  the  energy  transport  medium,  i.e.,  the  fluid  phase  flowing 
through  the  packed  bed.  Refrigerant- 12  (known  as  R-12,  or  dichloro-difluoro-methane) 
was  chosen  as  the  vapor  phase  because  it  is  a  highly  inert  and  stable  compound  whose 
critical  point  is  well  above  the  range  of  temperatures  considered  in  the  present  study. 
Mor  over,  it  has  a  relatively  high  vapor  density  making  it  capable  of  carrying  more  thermal 
energy,  and  thus  requiring  smaller  volume  flow  rate  for  a  certain  application,  than  the 
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Figure  2.1:  Schematic  diagram  of  the  problem 
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typical  gases  such  as  air.  Different  materials  will  be  considered  for  solid  particles  in  order 
to  determine  the  effect  of  different  thermophysical  parameters  on  the  associated  transport 
processes.  In  this  problem,  the  high  temperature  vapor  in  the  reservoir  is  suddenly  allowed 
to  flow  through  the  packed  bed  initially  filled  by  R  -12  which  is  also  initially  at  a  uniform 
temperature  and  pressure  throughout  the  packed  bed.  The  pressure  in  the  reservoir  is 
higher  than  the  initial  pressure  in  the  packed  bed  and  it  remains  unchanged.  A  uniform 
pressure  equal  to  the  initial  pressure  in  the  packed  bed  is  maintained  at  the  exit  of  the  bed  so 
that  there  is  always  a  pressure  difference  between  the  left  hand  and  the  right-hand  sides  of 
the  bed,  thus  sustaining  the  vapor  flow.  In  essence  this  is  a  forced  convection  flow  through 
a  packed  bed.  Therefore,  relatively  higher  particle  Reynolds  number  flows  than  the  Darcy 
regime  flows  are  considered.  Since  the  Darcy  formulation  would  fail  for  these  flows, 
Ergun-Forchheimer  relation  is  employed  instead,  in  order  to  account  for  the  inertial  effects 
as  well  as  the  viscous  effects.  The  problem  is  analyzed  for  the  case  in  which  the  two 
horizontal  walls  of  the  packed  bed  are  maintained  at  constant  temperature. 


2.2  Mathematical  model 

The  governing  equations  for  the  transport  processes  in  porous  media  are  developed 
here  by  using  the  "local  volume  averaging"  technique.  This  method  enables  us  to  express 
the  governing  equations  in  a  localized  macroscopic  sense  for  the  porous  media,  which 
would  otherwise  require  a  very  complex  microscopic  analysis  that  does  not  seem  to  be 
feasible  even  with  the  most  advanced  computational  facilities.  The  fundamentals  of  the 
application  of  this  technique  to  porous  media  are  outlined  by  Whitaker  (1977)  in  his  work 
on  drying  of  porous  media.  In  this  work  the  governing  conservation  equations  have  been 
derived  for  Darcy  flow  cases.  Gray  (1975)  also  presented  an  extension  on  the  application 


of  the  volume  averaging  technique  to  convective  and  diffusive  terms  in  the  governing 
equations.  The  volume  averaging  process  is  performed  by  associating  a  small  volume  V 
closed  by  a  surface  S  with  every  point  in  the  porous  medium.  The  two  important  averages 
of  a  quantity  frequently  encountered  in  the  governing  equations  are  the  so  called  local 
volume  average  and  the  intrinsic  phase  average.  The  local  volume  average  of  a  quantity  0 
associated  with  phase  ¥  is  defined  as: 

<  0  >  =  [  OdV  (2.D 
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while  the  intrinsic  phase  average  of  a  quantity  0  associated  with  phase  VF  is  defined  as: 
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where  Vq>  is  the  volume  associated  with  phase  M'.  With  this  formulation  it  becomes 
reasonable  to  represent  the  intrinsic  properties  such  as  the  temperature  and  the  density  of  a 
phase  by  the  intrinsic  phase  averaged  quantities,  and  to  represent  the  frequently  used 
"superficial  velocity"  by  the  local  volume  averaged  quantity. 


Considerable  information  regarding  the  application  of  the  spatial-averaging  theorem 
for  establishing  the  local  volume  averaged  form:*  of  the  conservation  equations  for  the 
convective  and  diffusive  flows  in  porous  media,  and  the  models  formed  by  using  this 
principle  for  such  problems  can  be  found  in  the  works  of  Vafai  and  Whitaker  (1986),  and 
Vafai  and  Tien  ( 1981).  The  same  methodology  is  used  in  the  present  work  for  establishing 
the  model  for  analyzing  the  energy  and  momentum  transfer  in  the  packed  bed  with  forced 
convective  flow.  The  major  assumptions  and  simplifications  that  are  employed  in  this 
study  are: 


(1)  Refrigerant- 12  used  as  the  vapor  phase  obeys  the  ideal  gas  equation  of  state,  and  the 


operation  range  is  chosen  such  that  it  is  always  in  the  superheated  form. 

(2)  Natural  convection  effects  are  negligible.  This  is  a  very  reasonable  assumption 
based  on  the  fact  that  the  bed  is  dominated  by  forced  convection. 

(3)  The  solid  phase  is  incompressible,  and  the  packed  bed  is  assumed  to  have  uniform 
porosity  and  to  be  isotropic. 

(4)  Boundary  and  variable  permeability  effects  are  neglected  in  the  vapor  phase 
momentum  equation. 

(5)  The  temperature  range  used  in  the  analysis  of  the  problem  is  relatively  small  (20K); 
therefore,  the  variation  of  physical  properties,  such  as  thermal  conductivity, 
viscosity  and  specific  heat  capacity,  with  temperature  is  neglected. 

(6)  Due  to  the  relatively  low  temperature  range  considered  in  the  present  study  the  inter- 
particle  and  intra-partide  radiation  heat  transfer  is  negligible. 

With  these  points  taken  into  consideration  the  volume-averaged  governing  equations 

can  be  established  in  the  following  form  by  applying  the  volume  averaging  technique  to  the 

fundamental  flow  and  energy  equations  in  a  porous  medium: 

Vapor  phase  continuity  equation: 

The  volume  averaged  form  of  the  vapor  continuity  equation  for  a  compressible  fluid 

takes  the  following  form: 

J“(e  <  P  v  >v)  +  V.(<  p  v  >v  <vv>)  =0 
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(2.3) 


Vapor  phase  momentum  equation: 


The  vapor  phase  momentum  equation  which  incorporates  the  inertia  effects  as  well  as 
viscous  effects  by  the  use  of  Ergun-Forchheimer  relation  can  be  written  as  (Vafai,  1984): 


V<Pv>v  =  - 


<  Pv  >v  Fe  _ 


k!/2 


,  -  -  .  Vv  \LV  _ 

[<vv  >.<  Vv>]|^r-— <vv> 


(2.4) 


In  this  equation,  the  first  term  on  the  right-hand  side  represents  the  inertia  effect;  or  the 
pressure  drop  due  to  form  drag  resistance  of  the  packed  bed,  while  the  second  term 
represents  the  pressure  drop  due  to  viscous  effects.  Without  the  first  term,  the  equation 
reduces  to  Darcy  equation. 


When  assumption  (2)  above  is  incorporated  into  this  equation,  the  vapor  flow 
essentially  becomes  one-dimensional  and  it  assumes  the  following  form: 


0<PV  >v 
3x 


<  pv  >v  Fe 


K 
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<uv  > 


(2.4a) 


Vapor  phase  energy  equation: 


d  <  Tv  > 

E<pv>vcPv ^ +  c p  <pv>v<vv>.V<Tv>1 


V.  { 1 k  ,„rr 7  <  T  ,  >v  }  +  h  „ ; a  „  ( <  T,  >s  -  <  T  „  >*) 


(2.5) 


This  equation  accounts  for  conduction  and  transient  storage  of  energy  in  the  vapor  phase  as 
well  as  the  transport  of  sensible  heat  by  the  vapor  motion.  Vapor-to-solid  heat  transfer  is 
modeled  by  a  convective  heat  transfer  term  by  the  use  of  a  fluid-to-solid  heat  transfer 
coefficient. 


Solid  phase  energy  equation: 


(l-e)pscp  =  V. {  k ^  v  <  T ,  >s }  -  h „ a„ ( <  T ,  >s  -  <  T  v  >«)  (2.6) 

This  equation  follows  similarly  as  in  the  case  of  the  vapor  phase  energy  equation  except 
that  there  is  no  sensible  energy  transport  due  to  motion. 


Vapor  phase  equation  of  state: 

From  the  ideal  gas  law,  the  equation  of  state  can  be  written  as: 
<PV  >v  =  <pv  >VRV  <TV  >v 


(2.7) 


In  the  equations  presented  above,  the  effective  thermal  conductivities  for  the  solid  and 
the  vapor  phases  were  modeled  in  the  following  form: 


^vcff  ^ 

^scff =  ( 1  ’£)  Ks 


(2.8) 


The  permeability  of  the  packed  bed,  Kv,  and  the  geometric  function  F  in  the  vapor 
momentum  equation  can  be  obtained  from  the  experimental  results  of  Ergun  (Ergun,  1952; 
Vafai,  1984)  in  terms  of  the  porosity,  e,  and  the  particle  diameter,  dp,  as  follows: 


Kv  = 


3 , 2 
£  dp 


1500  -  e) 


(2.9) 


F  = 


1.75 

VT50  e3/2 


(2.10) 


In  equations  (2.5)  and  (2.6),  aM  represents  the  specific  surface  area  of  the  packed  bed,  i.e., 
the  surface  area  of  the  solid  panicles  per  unit  bulk  volume  of  the  packed  bed.  Based  on 
geometrical  considerations,  this  ratio  can  be  expressed  in  terms  of  the  particle  diameter  and 
porosity  as  (Dullien,  1979): 

6(1  -  £) 

asv  ~  j 

(2.11) 

For  the  fluid-tc-particle  heat  transfer  coefficient,  hsv,  it  is  necessary  to  choose  an 
empirical  correlation  from  a  range  of  available  experimental  results.  Considerable  amount 
of  experimental  work  has  been  carried  out  for  determining  this  quantity  for  packed  beds  for 
different  sizes,  shapes  and  packing  configurations  of  solid  particles.  Typical  experimental 
investigations  and  some  reviews  on  fluid-to-particle  heat  transfer  coefficients  may  be  found 
in  the  works  of  Raumeister  and  Bennett  (1958),  Bhatincharyya  and  Pei  (1975)  and  Barker 
(1965).  Considering  the  ranges  of  the  panicle  diameter,  dp,  and  the  paiticle  Reynolds 
number,  Rep,  used  in  this  investigation,  the  empirical  correlations  established  by  Gamson 
et  al.  (1943)  were  found  to  be  appropriate  for  use  in  the  present  work.  These  empirical 
correlations  which  were  presented  in  terms  of  Colburn-Chilton  Jh  factors  t.an  be  expressed 

in  the  following  form  after  some  manipulation: 


hsv  =  1.  064  cpG 


f  dpGl 

l  k  1 

l  R  J 

-  0.  4! 


for 


dpG 

R 


>350  (turbulent)  (2.12) 


hsv  =  18. 


CpG 


fC|,H  "s 

rdPG) 

l  k  J 

l  R  J 

dpG 

R 


for 


<40  (laminar)  (2.13) 


where  G  represents  the  rate  of  mass  flow  through  a  unit  surface  area  perpendicular  to  the 
flow  direction. 


Thus  equations  (2.3)  through  (2.13)  form  the  model  for  analyzing  the  energy  and 
momentum  transport  in  unsteady  forced  convection  flow  of  a  gas  through  a  packed  bed. 
This  is  basically  a  system  of  five  coupled  governing  equations  which  must  be  solved  for 
analyzing  the  five  field  variables,  namely  <  pv  >v ,  <  u v  >,  <  Tv  >v ,  <  Ts  >s ,  and 
<  Pv  >v  • 


2.2.1  Initial  and  boundary  conditions 

In  the  problem  under  consideration,  the  initial  conditions  are  mathematically 
expressed  as: 

Pv(x,y,t=0)  =  P0 

Tv(x,y,t=0)  =  Ts(x,y,t=0)  =  T0  (2.14) 

uv(x,y,t=0)  =  0 

The  boundary  conditions  used  are  as  follows: 

Pv(x=0,y,i)  =  Pin,  Pv(x=L,y,t)  =  P0 

Tv(x=0,y,t)  =  Tin  att>0+  (2.15) 

Tv(x,y=0,t)  =  Tv(x,y=H,t)  =  Ts(x,y=0,t)  =  Ts(x,y=H,t)  =  T0 

Constant  wall  temperature  boundary  conditions  were  employed  in  order  to  study  the  two- 
dimensional  effects  in  the  transport  phenomena  in  the  packed  bed.  This  kind  of  boundary 
conditions  will  give  an  insight  about  the  behavior  of  the  packed  bed  during  the  energy 
release  cycle. 


To  analyze  the  effect  of  different  parameters  on  energy  and  momentum  transport  in 
the  packed  bed,  different  solid  particle  sizes  and  different  materials  for  the  solid  phase  were 
considered.  A  complete  set  of  computations  were  carried  out  for  particle  diameter  values  of 
1  mm,  2.5mm  and  5  mm  for  the  set  of  different  materials  which  were  considered  for  the 
solid  phase.  Several  runs  were  also  performed  for  a  particle  diameter  of  1  cm.  The  three 
different  materials  considered  for  the  solid  phase  in  this  study  were  lithium-nitrate- 
trihydrate,  sandstone  and  1%  Carbon-steel.  The  main  criterion  in  selecting  these  materials 
was  to  obtain  a  wide  range  of  thermophysical  properties  and  explore  their  effects  on  the 
qualitative  behavior  of  the  transport  phenomena  since  the  properties  of  the  bed  material  may 
differ  widely  from  one  application  to  another.  The  average  porosity  was  chosen  to  be  equal 
to  the  asymptotic  value  of  0.39  throughout  the  packed  bed.  This  is  a  valid  assumption  for 
all  cases  in  which  the  ratio  of  the  solid  particle  diameter  to  the  characteristic  global 
dimension  of  the  packed  bed  is  below  a  certain  value  as  determined  by  Benanati  and 
Brosilow  (1962).  The  following  physical  data  were  used  for  the  boundary  conditions  and 
for  the  properties  of  the  materials  considered  in  the  numerical  computations: 

Refriaerant-12  Lithium-nitrate-trihvdrate 

R  =  68.7588  J/kg.K  cp  =  2090  J/kg.K 

cp  =  602  J/kg.K  k  =  0.5  W/m.K 

k  =  0.0097  W/m.K  p  =  1 550  kg/m3 

p  =  12.6X1  O'6  Pa.s 


Sandstor 


1  %  Carbon-steel 


Cp  =  710  J/kg.K 
k  =  1.83  W/m.K 
p  =  2200  kg/m3 


cp  =  473  J/kg.K 
k  =  43  W/m.K 
p  =  7800  kg/m3 


P0=  100  kPa  Tin  =  300  K  To  =  280K 

The  effects  of  three  important  parameters  have  been  studied  extensively  in  this 
investigation.  These  are  the  particle  Reynolds  number,  the  Darcy  number  and  the  ratio  of 
the  solid  phase  thermal  diffusivity  to  liquid  phase  thermal  diffusivity.  Different  Reynolds 
numbers  were  obtained  by  varying  the  pressure  in  the  reservoir,  i.e.,  by  applying  different 
pressure  gradients  across  the  packed  bed.  Different  Darcy  numbers  resulted  from  different 
particle  sizes  since  the  permeability  of  the  packed  bed  varies  with  the  panicle  diameter  as 
given  in  equation  (2.9).  The  ratio  of  solid  to  vapor  thermal  diffusivity,  (ag/OLy),  was  of  a 
different  order  of  magnitude  for  each  of  the  three  materials  considered  for -the  solid  phase. 


2.3  The  numerical  code  and  solution  procedure 

A  numerical  code  was  developed  for  solving  the  system  of  governing  equations 
presented  in  the  previous  section.  The  Beam- Warming  type  implicit  method  which  is 
frequently  used  in  compressible  flow  problems  in  regular  media  was  not  applicable  in  this 
case  because  of  the  nature  of  the  equations.  Therefore,  it  was  found  necessary  to  employ 
explicit  schemes.  In  developing  the  finite  difference  forms  of  equations  (2.3)  through 
(2.7),  central  difference  approximations  were  used  for  most  of  the  spatial  derivatives  for 
the  inner  grid  points.  Spatial  derivatives  for  the  grid  points  on  the  left  and  right  boundaries 
were  formed  by  first-order  forward  and  backward  difference  approximations  respectively. 
Upwind  differencing  was  implemented  for  the  convective  terms  in  the  vapor  energy  and 
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continuity  equations  in  order  to  ensure  the  stability  of  the  numerical  scheme.  It  was  found 
through  numerical  experimentation  that  using  central  differencing  in  the  convective  term  of 
the  continuity  equation  did  not  make  any  appreciable  change  in  the  results. 

While  the  vapor  density  at  each  grid  point  on  the  left  boundary  was  determined  from 
the  equation  of  state  by  using  the  boundary  conditions  for  temperature  and  pressure,  the 
density  at  each  grid  point  on  the  right  boundary  was  determined  from  the  equation  of  state 
by  using  the  vapor  temperature  computed  from  equation  (2.5)  and  the  pressure  given  by  the 
right  boundary  condition.  For  the  rest  of  the  grid  points  the  density  was  computed  by  the 
vapor  continuity  equation.  The  velocity  <  uv  >  for  all  the  grid  points  except  those  on  the 

right  boundary  were  determined  by  the  vapor  phase  momentum  equation  which  is  basically 
a  quadratic  equation  in  <  u  v  > .  The  velocity  at  the  grid  points  on  the  right  boundary  was 

computed  by  linear  extrapolation  from  the  preceding  two  grid  points  in  x-direction. 
Temperatures  of  the  solid  and  fluid  phases  were  determined  by  the  associated  energy 
equations  while  the  pressure  was  determined  from  the  equation  of  state. 

The  numerical  computations  were  performed  on  a  CRAY  XMP/28.  As  a  result  of 
vectorization  the  CPU  time  required  on  the  CRAY  for  a  certain  run  using  a  grid  mesh  of  21 
x  11  was  approximately  30  times  smaller  than  the  CPU  time  required  for  the  same  run  on  a 
VAX  8550.  Due  to  the  nature  of  the  governing  equations,  small  time  steps  had  to  be 
employed  for  the  stability  of  the  numerical  scheme.  A  systematic  decrease  in  the  grid  size 
was  carried  out  and  the  corresponding  stable  At  was  employed  until  an  agreement  to 
within  1 9c  was  achieved  in  the  solutions  between  the  consecutive  grid  sizes,  and  the  effect 
of  reducing  the  At  further  did  not  make  any  appreciable  change  in  the  solutions. 

A  number  of  interesting  observations  were  made  on  the  stability  of  the  numerical 
code.  It  was  found  that  in  some  cases  the  time  step  used  in  the  early  stage  could  be 
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increased  in  the  later  stage  by  an  order  of  magnitude  without  any  influence  on  the  accuracy 
of  the  results.  It  was  observed  that  for  a  fixed  particle  Reynolds  number,  Rep,  decreasing 
the  Darcy  number,  Da,  consistently  increased  the  stability  of  the  numerical  scheme  in  both 
the  early  and  the  later  stages.  Decreasing  Da  at  a  fixed  Rep  requires  using  smaller  particle 
diameter  (hence  smaller  permeability)  and  a  larger  pressure  difference  across  the  packed 
bed. 

On  the  other  hand,  for  a  fixed  Da,  smaller  values  of  Rep  translate  into  smaller  mass 
flow  rates  and  smaller  amounts  of  thermal  energy  flowing  into  the  packed  bed  per  unit 
time.  Therefore,  for  small  values  of  Rep  the  CPU  time  required  for  simulating  the  complete 
charging  of  the  packed  bed  becomes  excessively  large.  For  instance,  for  the  case  when 
lithium-nitrate-trihydrate  is  the  solid  phase  and  Da=1.18xl0'8,  the  CRAY  CPU  time 
required  for  simulating  the  complete  charging  of  the  packed  bed  was  of  the  order  of  1  hour 
for  Rep=l  100.  For  this  specific  case  the  At  employed  in  the  computations  was  1.7647x10' 
5  sec.  for  the  early  stage  while  it  was  possible  to  push  it  up  to  1.7647xl04sec.  for  the  later 
stage.  The  CPU  time  was  estimated  to  increase  significantly  for  small  Reynolds  numbers. 
Therefore,  full  simulation  of  such  cases  is  very  costly  on  CRAY  and  would  be  extremely 
tedious  on  regular  main  frames.  For  this  reason  no  attempt  was  made  to  perform  the 
complete  simulation  for  such  cases. 


?. . 4  Benchmark  of  the  numerical  eode 

In  order  to  gain  confidence  in  the  accuracy  of  the  results  of  the  numerical  code 
developed,  it  was  necessary  to  perform  a  benchmark  against  some  analytical  solutions, 
experimental  results  or  other  reliable  numerical  results.  As  mentioned  earlier,  no  complete 
analysis  of  unsteady,  nonisothermal  and  nonequilibrium  flow  of  a  gas  through  a  porous 


medium  is  available  in  the  literature.  Therefore,  in  orde:  to  check  the  numerical  code  it  was 
necessary  to  compare  its  results  with  different  limiting  case  analytical  solutions  dealing  with 
different  aspects  of  energy  and  momentum  transport  in  porous  media.  This  benchmarking 
was  performed  in  two  main  sections.  These  were  namely  the  comparison  of  momentum 
transport  (pressure  and  hence  velocity  distributions)  in  an  isothermal  flow  of  an  ideal  gas 
through  a  porous  medium,  and  the  comparison  of  energy  transport  (temperature 
distribution  in  solid  and  fluid  phases)  of  an  incompressible  fluid  flow  through  a  porous 
medium  with  no  local  thermal  equilibrium  between  the  solid  and  vapor  phases.  The  two 
sources  chosen  for  this  benchmarking  were  the  most  relevant  sources  that  could  be  found 
for  possible  comparisons. 

The  benchmarking  of  the  energy  transport  was  carried  out  by  comparing  the  analytical 
solution  of  the  simplified  Schumann  model  presented  by  Riaz  (1977)  against  the  results 
generated  by  the  present  numerical  code.  This  was  achieved  by  making  the  necessary 
adjustments  to  our  model  in  older  to  reduce  it  to  a  system  which  was  equivalent  to  the 
simplified  Schumann  model,  i.e. ,  by  omitting  a  number  of  the  transport  terms.  This 
comparison  is  depicted  in  Figure  2.2  for  the  fluid  phase  and  solid  phase  temperature 
distributions  iespectivoly  in  terms  of  the  dimensionless  variables  that  appear  in  the  work  of 
Riaz  (1977).  As  may  be  seen  from  these  figures,  the  agreement  between  our  numerical 
results  and  the  analytical  solution  is  excellent. 

The  benchmarking  of  the  momentum  transport  was  carried  out  by  comparing  our 
results  with  the  analytical  solution  of  a  one-dimensional  isothermal  flow  of  an  ideal  gas 
through  a  semiinfinite  porous  medium  presented  by  Kidder  (1957).  Although  this 
analytical  solution  was  obtained  for  a  semiinfinite  porous  medium,  it  could  be  safely 
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Figure  2.2:  Comparison  of  the  temperature  solutions  of  the  present  work  against  the  analytical  solutions  presented  by  Riaz 
(1977) 


compared  with  our  numerical  results  for  small  times,  during  which  the  information  from 
the  boundary  had  not  reached  the  total  length  of  the  packed  bed  under  consideration.  The 
comparison  of  the  analytical  solution  with  the  results  generated  by  the  present  numerical 
code  is  shown  in  Figure  2.3  in  terms  of  the  dimensional  distance  and  time  and  the 
dimensionless  pressure  which  appears  in  the  work  of  Kidder  (1957).  Again,  a  very  good 
agreement  can  be  observed  between  our  numerical  results  and  the  analytical  solution. 


2.5  Results  and  discussions 

For  the  present  analysis,  it  was  found  to  be  appropriate  to  present  the  results  in  a  non- 
dimensionalized  form.  The  reference  quantities  that  are  used  in  nondimensionalizing  the 
variables  are  the  reference  density,  reference  velocity  and  reference  pressure  denoted  by 
p*,  u*  and  P*  respectively.  For  convenience,  the  reference  pressure,  P*,  for  each  run  is 
taken  to  be  the  average  of  the  initial  pressure  Pt)  and  the  pressure  of  the  vapor  flowing  from 
the  reservoir  into  the  packed  bed,  Pj„.  Likewise,  the  reference  density,  p*,  is  computed 
from  the  equation  of  state  by  using  P*  and  the  average  value  of  the  initial  and  the  vapor 
inflow  temperatures.  Finally  the  refeience  velocity,  u*,  is  computed  from  equation  (2.4a) 
by  using  p*  and  the  global  pressure  gradient  applied  across  the  packed  bed.  The 
dimensionless  time  and  tempei  ature  Jiosen  for  presenting  the  results  are  defined  as  t  - 
tu*/L,  and  0  =  (T-T0)/(Tin-T0)  respectively.  For  each  of  the  computational  runs  two 
distinct  regimes  were  observed.  These  two  regimes  will  be  referred  to  as  the  '"early 
stage"'  and  the  "'Liter  stage"’.  The  early  stage  lasts  for  a  very  short  period  of  time  which 


starts  with  the  application  of  the  high  temperature  and  high  pressure  boundary  conditions 
on  the  left,  and  ends  when  the  piessiire  distribution  across  the  packed  bed  becomes  almost 
linear.  This  process  takes  place  extremely  fast,  and  afterwards  the  pressure  distribution 
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Figure  2.3:  Comparison  of  the  pressure  solution  of  the  present  work  with  the  analytical  solution  of  Kidder  (1957) 


remains  in  a  quasi-steady  linear  state.  During  the  early  stage  no  appreciable  thermal 
penetration  depth  develops  in  the  packed  bed.  The  later  stage  is  the  period  of  time  following 
the  early  stage. 

The  vapor  velocity,  density  and  pressure  distributions  along  the  x-direction  at  the 
midplane  of  the  packed  bed  during  the  early  stage ,  for  the  case  where  the  solid  phase 
material  was  lithium-nitrate-trihydrate,  are  depicted  in  Figure  2.4.  It  was  observed  that  for 
this  case  the  variation  of  different  variables  in  the  y-direction  for  the  core  region  of  the 
packed  bed  did  not  exceed  1%.  Therefore,  it  is  very  reasonable  to  present  the  results  for 
this  case  in  one-dimensional  form.  The  variation  in  the  velocity  distribution  can  be 
explained  by  the  variation  of  the  pressure.  Initially  when  the  high  pressure  boundary 
condition  is  applied,  the  large  pressure  gradient  causes  high  velocities.  Gradually  as  the 
pressure  propagates  in  time  and  space,  the  large  pressure  gradients  die  away,  and  the 
velocities  follow  the  same  trend.  As  mentioned  earlier,  there  is  no  appreciable  thermal 
penetration  depth  during  the  early  stage.  Therefore,  there  is  a  sharp  decrease  in  the  vapor 
temperature  near  the  entrance  of  the  bed.  This  drop  in  temperature  is  more  pronounced 
than  the  corresponding  pressure  drop.  The  resulting  peak  value  of  the  vapor  density  close 
to  the  entrance  region,  therefore,  follows  from  the  equation  of  state.  Beyond  the  peak 
region  the  temperature  is  almost  constant  while  the  pressure  decreases  in  the  x-direction. 
Hence,  based  on  the  equation  of  state  the  vapor  density  will  follow  the  same  trend  as  the 
vapor  pressure.  On  the  other  hand,  the  temporal  increase  in  the  vapor  density  during  the 
early  stage  occurs  because  during  this  period  the  mass  flux  at  any  point  is  larger  than  the 
mass  flux  at  a  point  on  its  right,  i.e.,  there  is  a  continuous  net  mass  flux  into  the  bed. 
Therefore,  the  second  term  on  the  left-hand  side  of  equation  (2.3)  has  a  negative  value. 
Thus,  the  density  has  to  inctease  in  time  during  the  early  stage  in  order  to  satisfy  the 
continuity  equation  of  the  vapor  phase. 
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Figure  2  4:  Variations  of  different  field  variable  distributions  at  the  mid-plane  of  the 


packed  bed  during  the  early  stage 


The  variations  of  the  solid  and  vapor  temperatures,  and  the  density,  velocity  and 
pressure  of  the  vapor  phase  during  the  later  stage  for  the  above-mentioned  case  are  shown 
in  Figure  2.5.  It  should  be  noted  that  the  later  stage  starts  after  the  pressure  distribution  in 
the  packed  bed  becomes  almost  linear  and  remains  in  that  quasi-steady  state  thereafter. 
During  this  stage  the  thermal  penetration  depth  develops  as  shown  in  the  figure.  Since  the 
effect  of  the  transient  ierm  in  the  vapor  continuity  equation  becomes  much  less  pronounced 
at  the  end  of  the  early  stage,  an  almost  steady-state  continuity  equatio.i  prevails  in  the  later 
stage.  This  requires  the  existence  of  an  inverse  relationship  between  the  vapor  density  and 
vapor  velocity  distributions.  This  behavior  can  be  seen  in  the  figure  as  almost  a  mirror 
image  trend  in  the  variations  of  the  density  and  the  velocity,  i.e.,  as  one  increases  the  other 
one  decreases.  The  vapor  pressure  distribution  is  almost  linear  in  the  x-direction  and 
remains  unchanged.  Furthermore,  at  any  time  before  the  packed  bed  is  fully  charged,  the 
slope  of  the  temperature  profile  is  much  steeper  than  the  slope  of  the  pressure  profile 
within  a  major  portion  of  the  thermal  penetration  depth.  Therefore,  during  such  times 
‘he  vapor  density  increases  in  the  x-direction  within  this  zone.  Beyond  this  zone  the  vapor 
temperature  is  almost  constant  and  the  vapor  pressure  decreases  linearly,  and  so  the  vapor 
density  also  decreases  linearly.  The  advancement  of  the  thermal  penetration  front  which 
translates  into  a  temporal  increase  in  temperature  at  a  given  x-location  causes  a  decrease  in 
the  vapor  density  in  time  at  that  location.  This,  of  course,  follows  from  the  equation  of 
ctate.  When  ultimately  a  uniform  temperature  is  reached  throughout  the  packed  bed,  the 
densit)  distribution  becomes  similar  to  the  pressure  distribution  as  dictated  by  the  equation 
of  state.  From  the  variation  of  the  different  field  variables  in  Figure  2.5,  it  becomes  evident 
that  the  flow  field  is  very  much  influenced  by  the  temperature  field. 

A  typical  example  run  in  which  two-dimensional  variations  are  very  pronounced  is 
considered  next.  The  solid  particle  material  considered  for  this  case  was  1%  Carbon-steel. 
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Figure  2.5:  Variations  of  different  field  variable  distributions  at  the  mid-plane  of  the 
packed  bed  during  tire  later  stage 
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Figure  2.5  (continued) 
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Figures  2.6(a)  and  2.6(b)  depict  the  distributions  of  the  solid  and  vapor  phase 
temperatures,  the  vapor  density  and  the  vapor  pressure.  As  may  be  seen  from  these 
figures,  there  is  considerable  heat  loss  through  the  horizontal  walls  of  the  packed  bed  due 
to  conduction  heat  transfer  in  y-direction  and  a  strongly  two-dimensional  temperature 
distribution  exists  within  the  packed  bed.  Apparently,  a  one-dimensional  modeling  for 
such  a  case  will  yield  erroneous  results.  In  order  to  better  track  the  two-dimensional 
behavior  of  the  variables  a  grid  mesh  of  21  x  21  was  employed  for  the  cases  in  which 
sandstone  or  steel  was  considered  as  solid  material. 

It  was  observed  from  the  numerical  computations  that  during  the  later  stage  of  each 
run,  the  value  of  the  vapor-to- solid  heat  transfer  coefficient  computed  for  each  grid  point  at 
each  time  step  was  very  stable  and  did  not  differ  more  than  2%  from  the  reference 
vapor-to-solid  heat  transfer  coefficient  which  was  based  on  the  reference  density,  p*,  and 
the  reference  velocity,  u*.  This  suggests  that  an  assumption  that  takes  this  value  to  be  a 
constant  for  all  the  grid  points  during  the  later  stage  would  be  very  reasonable. 
Considerable  CPU  time  can  be  saved  in  this  way  by  eliminating  the  otherwise  required 
computations  from  eqns  (2.12)  and  (2.13)  at  each  time  step.  However,  all  of  the  runs  for 
this  investigation  were  performed  by  using  the  variable  heat  transfer  coefficient. 

2.5.1  Effect  of  the  Darcy  number  (Da); 

In  order  to  determine  the  effect  of  the  Darcy  number  on  energy  and  momentum 
transport  in  the  packed  bed,  different  cases  with  different  Darcy  numbers  were  compared. 
These  comparisons  were  performed  for  a  fixed  particle  Reynolds  number  of  1 100  using  the 
same  solid  material.  Again  lithium-nitrate-trihydrate  was  used  for  convenience  in 
presenting  the  results,  since  for  this  material  the  distribution  of  the  field  variables  were 
essentially  one-dimensional  in  all  of  these  cases.  Figure  2.7  depicts  the  temperature 
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Figure  2.7:  Temperature  distributions  in  the  vapor  and  solid  phases  for  different  Darcy 


numbers  for  a  fixed  particle  Reynolds  number 
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Figure  2.7  (continued) 


distributions  in  the  solid  and  vapor  phases  along  the  packed  bed  at  different  time  levels  for 
Darcy  numbers  of  1 .1 8xlO-8, 7.38xl0'8  and  2.95xl0‘7.  The  time  levels  in  this  figure  were 
chosen  such  that  they  would  result  in  a  reasonable  presentation  of  the  penetration  fronts  in 
each  case.  It  can  easily  be  seen  that  as  the  Darcy  number  decreases,  the  difference  between 
the  solid  and  vapor  phase  temperatures  decreases.  The  reason  for  this  is  that  a  decrease  in 
the  Darcy  number  translates  into  a  decrease  in  the  particle  diameter  (see  equation  (2.9)).  As 
the  particle  diameter  decreases,  the  specific  surface  area,  asv,  of  the  packed  bed  increases, 
thus,  increasing  the  vapor-to-solid  heat  transfer  interaction  by  offering  a  larger  surface  area. 
Hence,  the  exchange  of  heat  between  the  solid  and  vapor  phases  becomes  more  efficient. 
Therefore,  the  local  thermal  equilibrium  assumption  would  be  more  justifiable  at  low  Darcy 
number  flows. 

2.5.2  Effect  of  the  particle  Reynolds  number  (Repl; 

The  procedure  followed  for  determining  the  effect  of  Rep  on  transport  phenomena 
was  similar  to  the  one  carried  out  in  determining  the  effect  of  the  Darcy  number.  This  time 
Da  was  fixed  at  7,38.\10'8  for  the  same  solid  phase  material  while  Rep  was  varied.  The 
temperature  distributions  for  the  solid  and  vapor  phases  along  the  x-direction  at  the 
horizontal  midplane  of  the  packed  bed  are  shown  in  Figure  2.8  for  Rep  of  400,  1 100  and 
2126.  As  can  be  seen  from  these  figures,  the  difference  between  the  temperatures  of  the 
solid  and  fluid  phases  at  any  point  increases  with  an  increase  in  Rep.  This  increase  should 
be  attributed  to  the  veloeit)  of  the  How  since  the  specific  surface  area  of  the  packed  bed  is 
the  same  in  all  of  these  cases.  As  the  velocity  of  the  flow  increases,  i.e.,  Rep  increases,  the 
time  for  the  solid-to  vapor  heat  interaction  decreases.  This  will  cause  a  decrease  in  the 
efficiency  of  heat  c \c flange  between  the  solid  and  vapor  phases  and,  hence,  the  deviation 
from  the  local  thermal  equilibrium  will  increase. 
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Figure  2.8:  Temperature  distributions  in  the  vapor  and  solid  phases  for  different  particle 
Reynolds  numbers  for  a  fixed  Darcy  number 
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Figure  2.8  (continued) 
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The  combined  effect  of  increasing  both  Rep  and  Da  was  found  to  have,  as  expected,  a 
much  more  pronounced  effect  on  the  local  thermal  equilibrium  assumption.  It  becomes 
evident  that  for  high  Rep  and  Da  flows  the  local  thermal  equilibrium  assumption  becomes 
very  erroneous.  Similar  qualitative  behavior  was  also  found  to  be  valid  for  the  cases  in 
which  sandstone  or  steel  was  used  as  the  solid  phase.  Lithium-nitrate-trihydrate  was 
chosen  just  for  its  convenience  in  presenting  the  results. 

2.5.3  Effect  of.  the  thermal  diffusivities  ratio  (ah7ar); 

It  was  found  through  numerical  experimentations  that  the  ratio  of  the  solid  phase 
thermal  diffusivity  to  liquid  phase  thermal  diffusivity,  as/av,  was  the  most  appropriate 
parameter  in  representing  the  combined  effects  of  k,/kv,  cPs/cPv  and  ps/pv  ratios,  all  of 

which  are  determining  factors  in  the  overall  heat  transfer  process.  Furthermore,  from  the 
numerical  experimentations  it  was  also  found  that  one-dimensional  approach  was  very 
satisfactory  for  all  cases  in  which  lithium-nitrate-trihydrate  was  used  as  the  solid  phase  and 
for  which  aja.w  was  of  the  order  of  0.035  to  0.075,  which  resulted  from  different  vapor 
density  values  corresponding  to  the  average  operating  pressures  applied  for  different  runs. 
Therefore,  the  effect  of  o.Jav  will  be  shown  for  cases  in  which  sandstone  and  steel  were 
used  as  the  solid  phase.  The  tempeiature  distributions  in  the  solid  and  fluid  phases  for  two 
cases  in  which  sandstone  and  steel  were  used  respectively  are  shown  in  Figs.  2.9(a)  and 
2.9(b).  In  both  of  these  figures  the  value  of  Rep  was  1 100  and  the  Darcy  number  was 
chosen  to  be  1.1 8x1 0'8.  In  the  case  of  sandstone  as/av  was  equal  to  0.57  where  in  the 
case  of  steel  it  was  5.71.  As  can  be  seen  from  these  figures  the  temperature  distributions 
are  strongly  two  dimensional  in  the  case  of  steel  and  mildly  two-dimensional  in  the  case  of 
sandstone  for  which  one  dimensional  approach  may  gi\e  fairly  good  results.  It  should  be 
noted  that  the  same  qualitative  behavior  was  observed  for  all  ranges  of  Rep  and  Da 


Figure  2.9:  Temperature  distributions  in  the  solid  and  vapor  phases  at  the  mid-plane  of  the  packed  bed  for  Rep  =1100, 
Da  =  1.18xl0*8,  (a)  at  %  =  103.9  for  the  case  with  sandstone,  (b)  at  x  =  568.7  for  the  case  with  steel 


considered  in  this  investigation  which  proves  that  two-dimensionality  effecis  are  not 
sensitive  to  the  panicle  Reynolds  number  and  Darcy  number. 

For  a  qualitative  assessment  of  the  validity  of  local  thermal  equilibrium  assumption 
and  the  strength  of  the  two-dimensional  effects  in  the  energy  and  momentum  transport,  the 
findings  of  the  numerical  computations  will  be  presented  in  an  integrated  form.  Figure 
2.10  shows  this  assessment  for  the  local  thermal  equilibrium  assumption  for  the  case  of 
lithium-nitrate-trihydrate,  sandstone  and  steel  for  the  entire  ranges  of  Rep  and  Da  which 
were  considered  in  this  work.  It  should  be  noted  that  for  all  the  runs  there  is  a  difference 
between  the  solid  and  vapor  phase  temperatures  at  time  equal  to  zero,  and  therefore,  the 
present  assessment  is  based  on  the  later  times  during  which  the  thermal  penetration  depth 
is  developed  to  some  appreciable  position.  The  dividing  lines  in  this  figure  were  decided 
by  the  maximum  difference  between  the  solid  and  the  vapor  phase  temperatures  relative  to 
the  overall  temperature  range  (difference  between  the  inlet  vapor  temperature  and  the  initial 
temperature,  20  K).  For  the  qualitative  ratings,  the  percentage  difference  relative  to  the 
overall  temperature  range  falls  into  the  following  categories:  very  poor,  more  than  15%; 
poor,  10-15%;  fair,  5-10%;  good,  1-5%;  very  good,  less  than  1%.  From  Figure  2.10  it 
can  be  concluded  that  Da  is  the  most  influential  parameter  in  determining  the  validity  of 
local  thermal  equilibrium.  The  Reynolds  number  based  on  particle  diameter  is  also  a  very 
important  parameter  in  this  regard.  The  local  thermal  equilibrium  assumption  becomes 
more  viable  as  both  Rep  and  Da  decrease  for  the  reasons  which  were  explained  in  the 
sections  analyzing  the  effects  of  these  parameters.  This  is  clearly  seen  from  the  lower  left 
comer  of  each  figure.  By  the  same  token  it  is  reasonable  to  observe  the  opposite  behavior 
as  Rep  and  Da  increase,  i.e.,  approaching  the  top  right  corners  of  these  figures.  Also  it 
may  be  seen  that  the  therinoph>sical  parameters  have  a  much  less  pronounced  influence  on 
the  local  thermal  equilibrium  condition.  Only  for  the  case  of  steel  were  the  effects  of  the 
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Figure  2. 10:  Qualitative  assessment  of  the  local  thermal  equilibrium  for  the  cases  with 
lithium-nitrate  trihydrate,  sandstone  and  steel  as  solid  phase  respectively 
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Figure  2.10  (continued) 


physical  parameters  found  to  be  mildly  significant.  This  is  because  the  conduction  heat 
transfer  becomes  more  dominant  than  the  vapor-to-solid  convective  heat  transfer  in  this 
case.  Therefore,  the  conduction  heat  transfer  will  play  an  additional  role  in  the  case  of 
steel.  This  was  significantly  absent  in  the  cases  of  sandstone  and  lithium-nitrate-trihydrate. 

In  a  similar  manner  a  qualitative  assessment  of  the  strength  of  the  two-dimensional 
behavior  is  depicted  in  Figure  2.1 1.  This  figure  shows  the  behavior  for  a  fixed  Da  of 
7.38x10'^  for  the  entire  range  of  Rep  used  in  this  investigation,  and  the  behavior  for  a  fixed 
Rep  equal  to  400  for  the  entire  range  of  Da  used  in  this  investigation.  This  qualitative 
behavior,  for  each  case  in  the  figure,  was  found  to  be  exactly  the  same  for  the  entire  ranges 
of  the  corresponding  Da  or  Rep  which  were  considered  in  this  investigation.  Therefore,  for 
brevity  only  these  two  cases  are  presented.  The  above-mentioned  fact  and  Figure  2.11 
confirm  that  the  two-dimensionality  effects  are  not  sensitive  to  either  Rep  or  Da,  but  very 
sensitive  to  the  ajav  ratio,  i.e.,  to  the  thermophysical  properties  of  the  materials  which  are 
considered  for  the  bed. 


2.6  Conclusions 

Solution  of  the  transient  problem  with  step  change  boundary  conditions  showed  the 
existence  of  two  distinct  .stages,  namely  the  early  stage  and  the  later  stage.  It  was  found 
that  the  early  stage  is  mostly  dominated  by  changes  in  the  pressure  distribution,  and  the 
time  scale  for  the  changes  in  the  distribution  of  pressure  and  the  dependent  variables  is  very 
short  during  this  stage.  On  the  other  hand,  the  time  scale  for  changes  in  the  field  variables 
during  the  later  stage  was  found  to  be  relatively  larger,  and  that  development  of  the  thermal 
penetration  into  the  packed  bed  played  an  important  role  in  these  changes. 
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Figure  2. 1 1 :  Qualitative  assessment  of  the  two-dimensionality  effects 


Modeling  the  flow  of  the  vapor/gas  by  using  an  equation  of  state  was  found  to  be 
very  useful  since  this  way  it  was  possible  to  track  down  the  pressure  evolution  as  well  as 
the  variations  in  the  velocity  and  density  of  the  transporting  fluid  during  the  early  stage. 
These  can  be  crucial  information  during  the  start-up  procedures  of  certain  applications.  It 
was  also  shown  that  although  no  such  assumption  has  to  be  made,  the  pressure  distribution 
in  the  packed  bed  assumes  a  linear  form  during  the  later  stage. 

The  results  obtained  in  the  case  studies  performed  in  Section  II  were  also  helpful  in 
drawing  a  number  of  conclusions  on  the  qualitative  aspects  of  transport  phenomena  in 
packed  beds  generally  applicable  to  any  porous  medium.  The  first  significant  outcome  of 
these  analyses  was  related  to  the  local  thermal  equilibrium  (LTE)  assumption.  It  was 
clearly  seen  from  the  results  obtained  that  higher  particle  Reynolds  number  resulted  in 
larger  differences  between  the  solid  and  fluid  phases,  that  is  larger  deviations  from  local 
thermal  equilibrium.  Similar  behavior  was  seen  with  increase  in  the  Darcy  number  which 
translates  into  an  increase  in  the  particle  size  of  the  packed  bed.  On  the  other  hand,  the 
effect  of  the  thermophysical  properties  on  LTE  was  found  to  be  much  less  influential. 
Therefore,  it  can  be  concluded  that  the  validity  of  LTE  assumption  is  mostly  dependent  on 
the  speed  of  the  How  and  the  particle  size  of  the  porous  medium. 

The  two-dimensional  behavior  in  the  transport  processes  was  also  examined.  It  was 
found  that  while  the  two-dimensional  effects  are  very  sensitive  to  the  thermophysical 
properties  of  the  solid  and  fluid  phases  considered,  particle  size  or  speed  of  flow  do  not 
have  significant  effect  on  two-dimensional  behavior.  The  results  obtained  show  that  higher 
the  thermal  diffusivity  of  the  solid  phase  relative  to  that  of  the  working  fluid  the  more 
pronounced  are  the  two-dimensional  effects.  This  is  reasonable  since  with  higher  thermal 


diffusivity  in  the  solid  phase  the  diffusion  of  heat  in  the  transverse  direction  can  reach  the 
order  of  magnitude  of  the  convective  heat  transfer  between  the  solid  and  the  fluid  phases. 

Numerical  experimentations  showed  that  the  stability  of  the  numerical  code  for  the 
later  stage  simulations  was  primarily  dependent  on  the  magnitudes  of  the  particle  Reynolds 
number,  Rep,  and  Darcy  number,  Da.  The  maximum  stable  time  step,  At,  was  found  to 
increase  with  an  increase  in  Rep  or  decrease  in  Da.  While  the  reason  for  increased  stability 
with  decrease  in  Da  can  be  attributed  to  the  increased  damping  forces  in  the  packed  bed,  the 
reason  for  increase  in  the  stability  of  the  numerical  code  with  increase  in  Rep  should  be 
related  to  the  level  of  constraint  that  the  flowing  gas  is  under,  or  in  other  words  the 
magnitude  of  the  pressure  difference  applied  across  the  packed  bed. 

It  was  found  out  that  due  to  the  nature  of  the  governing  equations,  very  small  time 
steps  have  to  be  employed  in  the  numerical  simulation  of  the  transport  processes  in 
compressible  flow  through  porous  media.  Numerical  simulation  of  the  complete  thermal 
charging  process  of  the  packed  bed  for  cases  with  low  Rep  and/or  high  Da  requires 
excessive  CPU  time  even  on  vector  machines.  For  such  cases,  the  use  of  regular  main 
frames  becomes  highly  inefficient. 

One  significant  finding  in  Section  II  was  that  the  use  of  a  fluid-to-solid  heat  transfer 
coefficient  based  on  a  representative  reference  velocity,  i.e.,  a  nominal  Rep  and  a  reference 
density  throughout  the  later  stage  simulations  would  not  cause  a  significant  inaccuracy  in 
the  numerical  solutions.  This  is  of  engineering  importance  because  considerable  amount  of 
CPU  time  can  be  conserved  by  introducing  this  simplified  form  rather  than  performing 
computations  for  each  grid  point  at  each  time  step. 
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SECTION  III 


ANALYSIS  OF  MULTIPHASE  TRANSPORT  PHENOMENA 
IN  A  CONDENSING  FLOW  OF  A  VAPOR  IN  A  PACKED  BED 

In  this  section,  the  problem  analyzed  in  Section  II  is  extended  to  cover  the  multi¬ 
phase  transport  phenomena  in  a  packed  bed.  This  basically  involves  the  consideration  of 
the  same  flow  problem  at  a  pressure  range  in  which  the  vapor  (working  fluid)  density  may 
reach  the  saturation  vapor  density  causing  condensation.  A  rigorous  model  for  analyzing 
the  multiphase  transport  in  the  packed  bed  is  developed  by  considering  the  basic 
thermodynamics  of  the  condensation  piocess.  Numerical  solution  procedure  is  discussed 
in  detail.  This  is  followed  by  the  presentation  of  the  results  of  a  set  of  case  studies 
performed  for  this  problem. 


3 . 1  Statement  of  the  problem 

The  problem  under  consideration  is  the  transient  condensing  flow  of  a  vapor  through 
a  fixed  bed  of  regularly  sized  spheiical  solid  particles  packed  in  a  two-dimensional  channel. 
Figure  3. 1 .  depicts  the  schematic  diagram  of  the  problem  under  consideration.  The  packed 
bed  is  initially  filled  with  the  woiking  fluid  (Refrigerant- 12)  at  a  slightly  superheated  state 
and  the  whole  vapoi  and  solid  system  is  at  uniform  temperature  and  pressure.  R- 1 2  vapor 
at  a  higher  temperature  and  pressure  fiom  a  reservoir  is  suddenly  allowed  to  flow  through 
the  packed  bed,  thus  depositing  its  thermal  energy  to  the  solid  pai  tides  of  the  bed.  This 
basically  forms  a  step  change  in  the  temperature  and  pressure  at  the  inlet  boundary.  In  this 
problem  the  pressures  are  high  enough  so  that  at  certain  times  during  the  thermal  charging 


Figure  3.1:  Schematic  diagram  of  the  problem 
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of  the  packed  bed,  the  vapor  at  certain  locations  becomes  compressed  enough  so  that  its 
density  reaches  the  saturation  vapor  density  at  the  prevailing  temperature,  thus  causing 
condensation  to  take  place. 


3.2  Mathematical  model 

The  volume  averaging  technique  which  was  used  in  Section  II  will  be  employed  in 

this  section  too  in  developing  the  governing  equations  for  the  present  problem.  In 

establishing  the  model  for  the  analysis  of  this  problem,  the  following  assumptions  and 

simplifications  were  employed: 

(1)  The  width  of  the  packed  bed  is  significantly  larger  than  the  height  and,  therefore,  the 
problem  is  essentially  two-dimensional. 

(2)  The  vapor  phase  always  behaves  as  an  ideal  gas. 

(3)  There  is  no  local  thermal  equilibrium  (LTE)  between  the  solid  and  fluid  phases  but 
there  is  LTE  between  the  vapor  and  liquid  phases  when  there  is  condensate  present 
in  the  packed  bed. 

(4)  Nauru!  convection  effects  are  negligible  compared  to  the  forced  convection  in  the 
vapor  phase.  This  essentially  yields  a  one-dimensional  flow  in  the  vapor  phase. 

(5)  The  solid  and  the  liquid  phases  aie  incompressible,  and  the  packed  bed  has  uniform 
porosity  and  is  isotropic. 

(6)  Boundary  and  variable  porosity  effects  are  neglected. 

(7;  Inter-particle  and  intra-particle  radiation  heat  transfer  as  well  as  thermal  dispersion 
effects  are  neglected. 


(8)  Variations  of  physical  properties  such  as  thermal  conductivities,  specific  heat 

capacities,  viscosity  and  latent  heat  of  vaporization  with  temperature  are  neglected. 

With  all  these  assumptions  and  simplifications  taken  into  account,  the  volume 
averaged  governing  equations  can  be  established  in  the  following  form  (it  should  be  noted 
that  we  have  switched  to  subscripts  (3,  y,  a  and  f  which  denote  the  liquid,  vapor,  solid  and 
the  combined  fluid  (liquid+vapor)  phases  respectively): 

Vapor  phase  continuity  equation: 

y <  py  >Y)  +  V.(<  py  >Y<  \JY>)  =  -  <  fn  >  (3.1) 


Due  to  possible  phase  change  (condensation)  in  the  vapor  phase,  there  is  an  additional 
source  term,  <  m  > ,  in  this  equation  as  compared  to  equation  (2.3).  Moreover,  because 
of  phase  change,  the  volume  fractions  of  the  vapor  and  the  liquid  phases  will  be  functions 
of  time.  Therefore,  instead  of  a  constant  porosity  value,  e,  in  this  equation  we  have  a 
variable  for  the  vapor  phase  volume  fraction,  namely  e  Y(t) . 

Liquid  phase  continuity  equation: 


The  volume  averaged  continuity  equation  for  the  noncompressible  liquid  phase  can  be 
established  as: 
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(3.2) 


Vapor  phase  equation  of  motion: 


(3.3) 


<  p7  >' FeY 

v<py>  = - zm — [<V-<V]|^“3 
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k;<7v> 


This  equation  is  similar  to  the  one  used  in  Section  II  except  that  Ey  is  used  instead  of  £, 
and  the  definitions  of  F  and  Ky  are  slightly  different  as  will  be  discussed  later.  With 

assumption  (4)  above,  this  equation  reduces  to  the  following  one-dimensional  form: 


a  v  <pv>rFev  .  u  v 

£<Pi>= - —^2  ~  <ui'>2'“'k^<uv> 


(3.3a) 


Liquid  phase  equation  of  motion: 

The  liquid  phase  motion  is  formulated  by  Darcy’s  law  since  any  motion  of  the  liquid 
will  be  much  smaller  than  that  of  the  vapor  phase: 


k,PK 


<V  =  _X{k£VVk<T>V<Tf>t  +  (pp-<PY>I  >8} 


(3.4) 


Fluid  phase  energy  equation: 


r  t  -ia<Tr>r 

(C[i*Vcli'1  |i  ^  £7 <  Py  >  (ci))yj  ^  <m>Ah 
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+  [pp(Cp)  p  <  Vp  >  +  <  P7  >Y  (c p)  Y  <  Vy  >].  V  <  T f  >f  =  V.[k  fcff  v  <  T f  >f] 

+  h  Oft a  OP  [<  T  o  >0  - <  T  f  >f  ]  ■ +  h  ov 11  or  [<  T  O  >°  -  <  T  f  >f ] 


(3.5) 


Comparison  of  this  equation  with  equation  (2.5)  shows  the  additional  effect  of  the  liquid 
pha.,e  in  the  transient  stoiuge  teim,  conduction  term  and  the  sensible  heat  transport  term  due 
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to  motion.  Moreover,  there  is  an  additional  term  accounting  for  the  heat  transfer  due  to 
condensation  as  well  as  a  convection  heat  transfer  term  accounting  for  the  energy  transport 
between  the  solid  and  the  liquid  phases. 

Solid  phase  energy  equation: 

E0P<J(cp)l,~~5|2~  *  v-[k0crrv<T<j>°] 

“  ^1opaop[<  >  -<  >f]  ~  ^07  a  ot[<  ^  0  >  -<r^f>]  (3.6) 

In  this  equation  too,  the  energy  transport  term  between  the  liquid  and  the  solid  phases  is 
modeled  by  a  convective  heat  transfer  term. 

Volume  constraint  relation: 

Due  to  the  fact  that  the  liquid  and  vapor  volume  fractions  are  variable,  the  following 
volume  constraint  relation  ,  .ovides  a  coupling  equation,  namely: 

eo  +  eyOH  cp(t)  =  1  (3.7) 

Equation  of  state  for  vapor  phase: 

Following  assumption  (2),  the  equation  of  state  for  the  vapor  phase  becomes: 
<P7>Y=  <pY>YRT<Tf  >f  (3.8) 

Thermodynamic  relation  for  the  saturation  density  of  vapor: 
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exp 


(3.9) 


where  A  and  B  are  constants,  Tf  is  in  degrees  Kelvin  and  P  y,  s  is  in  kg/m3.  This  relation 
was  obtained  by  a  least  squares  curve  fitting  method  to  a  range  of  data  relevant  to  the 
present  study. 

Equations  (3.1)  through  (3.9)  yield  nine  equations  in  nine  unknowns,  namely  e^t) , 
e  y(  t ) ,  <p7>y,  <V,  <v|5>,  <PY>7,  <Tf>f,  <  T0  >° ,  and  <  m  > . 


For  the  present  analysis  the  effective  theimal  conductivities  for  the  solid  and  fluid 
phases  were  modeled  as: 

k(Teff  “  CakO 

(3.10) 

kfeff  =  ^y  +  epkp 

The  permeability  of  the  packed  bed  of  spherical  particles  is  given  in  the  following 
form  (Ergun,  1952): 


K  = 


e3dp 


150(1  -e)' 


(3.11) 


where  e  is  the  porosit)  and  d  is  the  panicle  diameter.  Therefoie,  the  permeability  for  the 
vapor  phase,  Ky,  and  the  geometric  factoi,  F,  in  the  vapor  phase  momentum  equation  can 
be  expressed  as  functions  of  dp  and  c.,  as  (Vafai,  1984): 


Ky  = 


3  .2 

£  y  Cl  p 


150(1  ~EyY 
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(3.12) 


For  the  relative  permeability  of  the  liquid  phase,  the  one  suggested  by  Udell  and  Fitch 


(1985)  was  used.  This  is  modeled  in  the  following  form: 


where: 


S  = 


s-s. 
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1  -S. 
■m 


and 


(3.14) 


where  S  is  the  normalized  saturation,  s  is  the  absolute  saturation,  and  sjm  is  the  "immobile" 
saturation.  The  immobile  saturation,  sjm,  is  the  critical  value  of  the  absolute  saturation.  If 
s  is  greater  •.  -i  oim  en  the  liquid  phase  in  the  packed  bed  becomes  mobile.  For  values  of 
s  lower  than  s|m,  any  liquid  phase  present  in  the  porous  system  will  be  immobile,  or  in  the 
so  called  pendular  state.  The  value  of  0. 1  used  for  sim  in  the  work  of  Kaviany  and  Mittal 

(1987)  will  be  used  in  the  present  work  because  of  lack  of  any  better  experimental  finding. 
With  this  value  of  sim,  the  critical  value  of  the  liquid  fraction  at  which  liquid  phase  becomes 
mobile,  i.e.,  £pcrjl,  was  computed  from  equation  (3.14)  to  be  0.039,  with  the  value  of 
porosity  of  the  packed  bed  taken  to  be  0.39,  which  is  the  average  asymptotic  value  for 
packed  beds  in  which  the  particle  diameter  to  packed  bed  diameter  is  below  a  certain  value 
(Benanati  and  Brosilow,  1962).  For  all  the  cases  studied  in  the  present  investigation  the 
maximum  value  of  Cp  never  reached  £(5Cril  and,  therefore,  the  liquid  was  always  immobile. 

It  should,  however,  be  noted  that  the  liquid  was  not  assumed  to  be  immobile  in  modeling 
the  problem,  but  the  fact  that  it  turned  out  to  be  immobile  was  the  consequence  of  the 
governing  physical  conditions  of  the  problem. 


The  empirical  correlations  for  the  fiuid-to-solid  heat  transfer  coefficients  that  were 


used  in  this  section  are  the  same  as  those  in  Section  II,  namely: 


haj  =  1.064  (cp)j 


rcpu- 
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hoj  =  18.1(cp)_ 
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<  40  (laminar) 


where  G  represents  the  mass  flow  rate  through  a  unit  surface  area  perpendicular  to  the 
direction  of  flow,  and  the  subscript  j  denotes  (5  or  y  for  the  liquid  or  vapor  phase 
respectively. 

The  specific  surface  area  of  the  packed  bed  for  the  vapor  phase  is  expressed  in  the 
following  form: 


6(1  —  eY  —  e  ) 


a  cry 


(3.16) 


Strictly  speaking,  this  cori elation  was  originally  derived  from  geometric  considerations  for 
a  fully  saturated  packed  bed  of  spherical  particles  for  a  single  fluid  phase,  in  the  form 
a  =6(1  -  £)  /dp  (Dull ien,  1979).  However,  since  £p  is  very  small  (less  than  0.01) 
computed  to  (approximately  0.38-0.39),  equation  (3.16)  provides  a  very  accurate 
estimate  for  a^ 


Also  from  an  analysis  of  the  representative  length  scales  and  volume  scales  of  the 
liquid  and  vapor  phases  of  the  working  fluid,  one  may  obtain  a  relation  between  the 
specific  surface  areas  a^y  and  aap  as: 
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and  hence: 
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Since  the  liquid  phase  was  always  in  pendular  state  and,  therefore,  hap  was  always 
zero,  due  to  the  fact  that  for  the  liquid  phase  the  mass  flow  rate  was  zero,  there  was  no  real 
need  for  the  above  information.  However,  aap  was  modeled  as  given  in  equation  (3.18) 
for  the  sake  of  complete  modeling. 


3.2.1  Initial  and  boundary  conditions 

The  physical  conditions  initially  prevailing  in  the  packed  bed  are  the  same  as  those  of 
the  previous  section.  They  are,  therefore,  mathematically  expressed  as  follows: 

Tj(x,y,t=0)  =  Ta(x,y,t=0)  =  T0 

PY(x,y,t=0)  =  P0  (3.19) 

Uy(x,y,t=0)  =  0 

The  pressure  on  the  right  boundary  is  kept  at  a  value  equal  to  the  initial  pressure  in  the 
packed  bed  while  vapor  at  a  fixed  high  temperature  and  pressure  is  supplied  at  the  inlet  of 


the  packed  bed.  The  mathematical  form  of  the  left  and  right  boundary  conditions  is  as 
follows: 


Tf(x=0,y,t)  =  Tin 

Py(x=0,y,i)  =  Pin  at  t>0+  (3.20) 

P7(x=L,y,t)  =  Pout  =  P0 

In  the  present  problem,  both  insulated  top  and  bottom  wall  boundary  conditions  as 
well  as  constant  temperature  top  and  bottom  wall  boundary  conditions  will  be  investigated. 
For  the  case  with  insulated  boundary  conditions  we  have: 


fcff 


dy 


=  0 


and 
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y  =  0,  y  =  H 


'acff  Sy” 
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y  =  0,  y  =  H 


(3.21) 


For  the  case  with  constant  temperature  boundary  conditions  at  the  top  and  bottom  walls,  the 
boundary  conditions  can  be  expressed  as: 


Tf(x,y=0,t)  =  T0(x,y=O,t)  =  Tbo[  =  TQ 
Tf(x,y=H,t)  =  Tc(x,y=H,t)  =  T,op  =  T0 


(3.22) 


3.2.2  Physical  conditions  for  the  numerical  runs 

In  the  present  investigation  too,  in  order  to  analyze  the  effects  of  different 
characteristic  parameters,  such  as  Rep  and  Da,  different  particle  diameters  ranging  from  1 
mm  to  5  mm,,  and  different  ranges  of  operating  pressures  were  considered.  Also  two 
different  solid  materials  were  considered  in  the  present  problem  for  determining  the  effect 
of  the  thermal  properties  of  the  solid  phase  on  transport  phenomena.  These  materials  were 


lead  and  1%  Carbon-steel.  The  value  of  the  average  porosity  of  the  packed  bed  was  taken 
to  be  0.39  following  the  arguments  in  Section  II.  The  numerical  values  of  different 
parameters  used  as  initial  and  boundary  conditions  are  given  as  follows: 

T0  =  300  P0  =  796  kPa  T^  =  350  K  Pout  =  796  kPa 

and  the  value  of  Pin  was  different  for  different  cases  ranging  from  800  kPa  to  866  kPa. 
The  numerical  values  of  the  other  physical  data  used  in  the  numerical  computations  are  as 
follows: 


R-12: 


ky=0.0097W/m.K  (cp)y=7 10.  J/kg.K  py=12.6xl0-6kg/m.s  Ry=0.0687588  J/kg.K 
Ahvap=111.3xl03J/kg  kp=0.0545  W/m.K  (cp)p=1115J/kg.K  ^=179.2xl0‘6kg/m.s 

pp=l  190.35  kg/m3 
Lead: 

kc  =  35.  W/m.K  pa  =  11340.  kg/m3  (cp)a  =  129.  J/kg.K 
1  %  Carbon  steel: 


kc  =  43.  W/m.K  p0  =  7800.  kg/m3  (cp)a  =  473.  J/kg.K 

3.3  Solution  procedure 

The  finite  difference  techniques  utilized  in  the  solution  of  the  problem  in  Section  II 
were  also  employed  in  the  solution  of  the  present  problem.  In  this  problem,  however, 
depending  on  whether  phase  change  (condensation)  takes  place  in  the  vapor  phase  or  not, 
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the  set  of  unknown  variables  and  governing  equations  changes.  The  vapor  phase 
continuity  equation  plays  an  important  role  in  determining  the  solution  format  of  the 
governing  equations.  The  basic  criterion  which  governs  the  phase  change  (condensation) 
in  the  working  fluid  is  the  attainment  of  the  saturation  vapor  density.  At  any  point  where 
the  vapor  density  reaches  the  saturation  density  corresponding  to  the  temperature  at  that 
point,  condensation  will  occur.  The  solution  format  of  the  condensation  mode  includes  a 
switching  from  the  vapor  phase  continuity  equation  to  the  thermodynamic  relation  giving 
the  saturated  vapor  density  for  the  solution  of  the  vapor  density,  while  the  vapor  continuity 
equation  is  then  used  for  the  solution  of  the  condensation  rate. 

The  stability  of  the  numerical  scheme  was  ensured  by  choosing  a  proper  combination 
of  Ax,  Ay  and  At.  A  systematic  decrease  in  the  grid  size  was  employed  for  obtaining  the 
convergence  of  the  numerical  scheme,  and  the  corresponding  stable  At  was  employed.  A 
compromise,  however,  had  to  be  made  between  the  accuracy  and  the  computer  CPU  time 
required  for  the  computational  runs.  A  41  x  21  grid  configuration  (which  gives  a 
dimensionless  Ax  of  0.025)  was  found  to  yield  qualitatively  and  quantitatively  good  results 
for  the  condensation  period  and  very  good  results  for  the  problem  for  the  cases  with 
insulated  boundary'  conditions.  A  41  x  41  grid  configuration  was  employed  for  the  cases 
with  constant  wall  temperature  boundaries. 


3.4  Results  and  discussions 

In  this  section,  part  of  the  results  that  will  be  presented  for  the  computational  runs 

Y  Y 

performed  will  be  in  terms  of  nondimensionalized  variables.  Variables  <  Py  >  ,  <Py> 
and  <  uY>  are  nondimensionalized  by  using  the  corresponding  reference  quantities, 
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namely  P*,  p*,  and  v*.  P*  was  chosen  to  be  100  kPa,  p*  was  calculated  from  the 
equation  of  state  by  using  P*  and  the  initial  temperature,  T0.  The  reference  velocity,  v*, 

was  computed  from  the  vapor  phase  momentum  equation,  (3.3a),  using  a  pressure  gradient 
which  was  based  on  the  global  pressure  difference  applied  across  the  packed  bed  and  a 
density  which  was  calculated  from  the  equation  of  state  by  using  To  and  Pav  (the  average  of 

the  inlet  and  exit  pressures).  Temperatures  of  the  solid  and  fluid  phases  are 
nondimensionalized  in  the  form  0=  (T-T0)  /(T^  -T0).  Time,  t,  is  kept  in 

dimensional  form  for  giving  an  insight  of  the  actual  magnitudes  of  the  durations  involved. 
Same  thing  is  done  for  the  condensation  rate  data  and  the  total  condensate  variations  as  well 
as  for  the  thermal  charging  data  of  the  packed  bed. 

As  in  the  case  of  the  solution  of  the  problem  in  Section  II,  two  distinctly  noticeable 

stages  were  observed  in  the  solution  of  the  problem  considered.  These  were  namely  the 

early  stage  and  the  later  stage.  The  early  stage  usually  lasts  for  a  very  short  period  of  time 

during  which  sharp  changes  in  the  distribution  of  certain  field  variables,  such  as  <  Py  >  , 
Y 

<  Py  >  and  ,  <  u y  > ,  occur  due  to  the  step  change  boundary  conditions  which  cause  very 
strong  transient  behavior. 

3.4.1  Problem  with  insulated  wall  boundary  conditions 

Figure  3.2  depicts  the  variations  of  the  density,  velocity  and  pressure  of  the  vapor  as 
well  as  the  variation  of  the  liquid  fraction  during  the  early  stage  at  the  midplane  of  the 
packed  bed  for  insulated  wall  conditions  for  a  case  in  which  lead  was  used  as  the  solid 
phase.  It  was  found  that,  for  cases  with  insulated  boundary  conditions,  one-dimensional 
formulation  would  be  very  accurate  since  there  was  no  appreciable  variation  of  the  field 
variables  in  the  y-direction.  Figure  3.2(a)  shows  that  the  vapor  pressure  distribution 
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Figure  3.2:  Variations  of  different  field  variable  distributions  during  the  early  stage 


evolves  and  becomes  almost  linear  during  the  early  stage.  There  is  no  appreciable  thermal 
penetration  during  this  period  and  the  spatial  variation  of  the  vapor  density  follows  the 
same  trend  as  the  vapor  pressure  as  dictated  by  the  equation  of  state,  except  at  the  entrance 
of  the  packed  bed.  The  tempora'  increase  in  the  density  is  a  result  of  the  transient  effects 
dictated  by  the  vapor  phase  continuity  equation.  At  points  where  the  vapor  density  reaches 
the  saturation  vapor  density,  condensation  occurs  and  the  liquid  phase  accumulates.  Due  to 
the  physical  conditions  of  the  problem,  the  period  of  time  during  which  the  vapor  phase  in 
the  packed  bed  maintains  a  high  enough  density  for  condensation  to  occur  is  relatively  short 
compared  to  the  total  thermal  charging  period  of  the  packed  bed.  The  variation  of  the  liquid 
fraction  is  shown  in  Figure  3.2(d).  In  Figure  3.2  the  early  stage  was  somewhat  extended 
to  include  the  time  during  which  more  than  99%  of  the  condensation  took  place  for  this 
case. 


Beyond  the  early  stage,  the  changes  in  the  field  variables  are  mainly  caused  by  the 
development  of  the  thermal  penetration  depth  since  the  pressure  distribution  remains  almost 
unchanged.  Distributions  of  the  field  variables  of  interest  by  time  during  the  later  stage  are 
shown  in  Figure  3.3.  In  Figure  3.3(a)  the  solid  lines  depict  the  solid  phase  temperature 
distribution  while  the  dotted  lines  depict  the  fluid  phase  temperature  distribution.  During 
the  early  stage,  the  effect  of  the  transient  term  as  well  as  the  condensation  (source)  term  in 
the  vapor  continuity  equation  die  out  and,  therefore,  during  the  later  stage  the  vapor  density 
variation  in  time  becomes  dependent  on  the  convective  term  in  this  equation.  The  mass 
flow  rate  in  the  packed  bed  becomes  constant  requiring  an  inverse  relationship  between  the 
vapor  density  and  velocity  Figures  3.3(c,d)  clearly  show  this  behavior  as  a  mirror  image 
type  of  tieiui  in  the  variations  of  the  vapor  density  and  vapor  velocity  at  any  instant.  The 
spatial  variation  in  the  \  apor  density,  on  the  other  hand,  can  be  explained  from  the  equation 
of  state.  At  any  instant  before  the  packed  bed  is  fully  charged,  the  slope  of  the  vapor 


temperature  profile  will  be  larger  than  the  slope  of  the  vapor  pressure  profile  for  a  certain 
length  of  the  packed  bed  (meaning  a  sharper  decrease  in  temperature  than  in  pressure) 
causing  an  increase  in  the  vapor  density,  where  exactly  the  opposite  behavior  can  be  seen  at 
the  locations  where  the  slope  of  the  vapor  temperature  profile  is  smaller  than  that  of  the 
vapor  pressure  profile  (meaning  a  sharper  decrease  in  pressure  than  in  temperature).  When 
the  packed  bed  becomes  thermally  fully  charged,  the  vapor  density  variation  follows 
exactly  the  same  trend  as  that  of  the  vapor  pressure  as  dictated  by  the  equation  of  state. 

The  overall  condensation  rate  in  the  packed  bed  was  computed  by  integrating  the 
individual  condensation  rates  at  all  the  grid  points  over  the  associated  volumes  at  each  time 
step.  The  average  overall  condensation  rate  per  unit  width  of  packed  bed  for  the  case  for 
which  the  early  and  later  stage  results  have  been  presented  is  depicted  in  Figure  3.4(a).  It 
can  be  seen  that  the  overall  condensation  rate  is  higher  at  the  beginning  when  the  transient 
effects  are  very  strong,  and  dies  away  as  the  sharp  changes  in  the  vapor  density  variation 
erode.  The  accumulative  condensate  in  the  packed  bed  was  also  computed  by  integrating 
the  condensation  rates  at  all  grid  points  over  the  associated  volumes  at  each  time  step  and 
totaling  with  the  previous  sum.  Figure  3.4(b)  shows  the  variation  in  time  of  the  amount  of 
total  condensate  in  the  packed  bed  per  unit  width  of  the  bed.  As  may  be  seen  from  this 
figure,  the  accumulation  is  fast  at  the  start  due  to  high  condensation  rate,  and  builds  up 
quickly  reaching  an  asymptotic  value  in  a  short  time. 

The  next  important  analysis  carried  out  was  related  to  the  thermal  charging  of  the 
packed  bed.  For  this,  we  needed  to  track  down  the  variations  of  the  relevant  parameters 
including  heat  flow  rates  and  energy  storage  rate.  The  variation  in  the  amount  of  thermal 
energy  flowing  into  and  out  of  the  packed  bed  per  unit  width  is  depicted  as  a  function  of 
time  in  Figure  3.5(a).  These  were  computed  by  integrating  the  energy  fluxes  of  the  vapor 
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Figure  3.4:  (a)  Variation  of  average  overall  condensation  rate  in  the  packed  bed 


(b)  Variation  of  total  accumulative  condensate  in  the  packed  bed 


(a) 


Figure  3.5:  (a)  Rate  of  heat  flow  into  and  out  of  the  packed  bed 


(b)  Thermal  charging  of  the  packed  bed  (insulated  top  and  bottom  walls) 


flowing  into  and  out  of  the  packed  bed  over  the  inlet  and  the  exit  crosssections.  A  very 
short  section  of  the  time  history  of  the  thermal  charging  process  at  the  beginning  was  left 
out  in  this  figure  in  order  to  obtain  a  better  scale  that  shows  the  variation  of  the  heat  flow 
rates  clearly  for  the  whole  charging  duration.  During  this  very  short  time  interval  the  rate 
of  heat  flow  into  the  packed  bed  was  several  times  larger  due  to  the  presence  of  large 
velocities  at  the  onset  of  the  charging  process,  i.e.,  during  the  early  stage.  Figure  3.5(b) 
shows  the  net  energy  stored  per  unit  width  of  the  packed  bed  as  a  function  of  time.  This 
was  computed  in  the  numerical  code  by  integrating  the  net  heat  flow  rate  into  the  packed 
bed  over  time.  The  value  of  the  total  net  energy  stored  in  the  packed  bed  was  also 
determined  by  a  thermodynamic  balance  analysis  between  the  initial  and  the  final  (fully 
charged)  states  of  the  packed  bed.  The  result  of  the  numerical  computation  was  found  to  be 
in  veiy  good  agreement  with  this  analytical  result. 


3.4.2  Problem  with  constant  temperature  wall  boundary  conditions 

The  same  problem  was  solved  for  constant  temperature  top  and  bottom  wall  boundary 
conditions.  As  expected,  strong  two-dimensional  behavior  was  found  in  the  variations  of 
many  of  the  field  variables.  Therefore,  rather  than  presenting  the  results  in  the  format  of 
the  previous  section  we  w  ill  resort  to  two-dimensional  contour  plots  for  the  field  variable 
distributions.  It  should  be  noted  that  the  solution  of  this  problem  has  early  stage  and  later 
stage  parts  just  as  in  the  case  of  insulated  boundary  conditions  case.  However,  the 
distribution  of  three  field  variables  at  three  time  levels  in  the  later  stage  will  be  presented 
here  for  convenience.  These  will  be  sufficient  to  show  the  two-dimensional  behavior  of  the 
problem  which  becomes  apparent  during  the  later  stage. 


Figure  3.6  depicts  the  distributions  of  the  fluid  temperature,  solid  temperature  and  the 
vapor  density  in  the  packed  bed  at  three  different  time  levels.  It  can  be  seen  from  this 
figure  that,  although  the  two-dimensional  behavior  of  the  field  variables  is  not  very 
significant  at  the  beginning,  it  becomes  very  much  pronounced  as  the  thermal  penetration  in 
the  bed  advances.  In  the  core  region  of  the  packed  bed  the  advancement  of  the  thermal 
front  follows  a  similar  trend  as  in  the  case  of  insulated  boundary  conditions.  Near  the  top 
and  bottom  walls,  however,  there  is  a  temperature  gradient  in  the  y-direction  due  to  heat 
loss.  Since  there  is  no  significant  variation  in  the  vapor  pressure  in  the  y-direction,  the 
vapor  density  variation  in  this  direction  is  primarily  determined  by  the  fluid  temperature 
variation.  Hence,  at  locations  closer  to  the  top  and  bottom  walls  where  the  fluid 
temperature  becomes  lower  the  vapor  density  becomes  higher.  In  the  core  region  of  the 
packed  bed  the  variation  of  the  vapor  density  in  x-direction  depends  on  the  slopes  of  the 
temperature  and  pressure  distributions  in  this  direction.  At  points  where  the  slope  of  the 
temperature  distribution  is  sharper  than  the  slope  of  the  pressure  distribution,  the  vapor 
density  will  be  increasing  and  vice  versa. 

The  average  overall  condensation  rate  per  unit  width  of  the  packed  bed  was  computed 
in  the  same  manner  described  in  the  previous  section.  Figure  3.7(a)  shows  the  variation  of 
this  quantity  in  time.  Upon  comparing  Figure  3.7(a)  with  Figure  3.7(b),  one  can  see  that 
condensation  lasts  slightly  longer  in  the  case  of  constant  wall  temperature  case.  This  is 
reasonable  because  due  to  the  heat  loss  from  the  top  and  bottom,  for  the  points  next  to  the 
top  and  bottom  boundaries  near  the  inlet  of  the  packed  bed,  it  takes  longer  for  the  vapor 
phase  to  reach  a  temperature  at  which  the  saturation  vapor  density  exceeds  the  vapor 
density.  Hence,  the  duration  of  condensation  at  these  points  is  prolonged.  This  also  gives 
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3.6:  Distribution  of  field  variables  in  the  packed  bed  at  different  time  levels  during  the  later  stage 
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Figure  3.7:  (a)  Variation  of  average  overall  condensation  rate  in  the  packed  bed 


(b)  Variation  of  total  accumulate  condensate  in  the  packed  bed 


rise  to  a  slightly  larger  amount  of  total  condensate  in  the  packed  bed.  This  can  be  seen  by 
comparing  Figure  3.7(b)  with  Figure  3.4(b). 

The  variations  in  the  amounts  of  heat  flowing  into  and  out  of  the  packed  bed  per  unit 
width  at  the  inlet  and  exit  are  shown  in  Figure  3.8(a)  as  functions  of  time.  In  this  case  too, 
a  small  time  slice  at  the  beginning  of  the  charging  process  is  omitted  from  this  figure  in 
order  to  obtain  a  better  scale  in  the  figure.  Figure  3.8(b)  depicts  the  energy  input  into  the 
packed  bed,  the  energy  that  is  lost  from  the  top  and  bottom  walls  by  conduction,  and  the 
net  energy  stored  in  the  packed  bed  per  unit  width  of  the  bed. 

3.4.2. 1  Effect  of  particle  Reynolds  number  (Rep)  on  condensation 

The  effect  of  Rep  on  condensation  was  investigated  by  running  three  cases  with 
different  Rep's  while  Darcy  number  (Da)  was  kept  constant  by  keeping  the  particle  diameter 
fixed.  The  solid  phase  in  all  these  cases  was  lead.  Different  Rep  values  were  obtained  by 
applying  different  inlet  pressure  boundary  condition  for  each  case.  Figure  3.9(a)  depicts 
the  variation  of  the  average  overall  condensation  rate  per  unit  width  of  the  packed  bed  for 
the  cases  in  which  the  nominal  Rep  values  were  500,  1,000  and  1,500  respectively.  The 
corresponding  variations  of  the  total  condensate  in  the  packed  bed  per  unit  width  as 
functions  of  time  are  shown  in  Figure  3.9(b).  From  these  figures  it  becomes  apparent  that 
the  higher  the  Rep  the  higher  will  be  the  condensation  rate  and  the  total  condensate 
accumulation.  On  the  other  hand,  higher  Rep  will  cause  slightly  shorter  duration  of 
condensation.  Higher  Rep  indicates  higher  mass  flow  rates  as  well  as  higher  pressure 
difference  applied  across  the  packed  bed  (larger  compression  forces  in  the  vapor  phase), 
thus,  higher  condensation  rate.  On  the  other  hand,  higher  Rep,  due  to  higher  pressure 
gradient  applied  means  potential  for  faster  propagation  of  the  pressure,  density  and 
temperature  distributions  in  the  vapor  phase  and,  thus,  shorter  condensation  time. 
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Figure  3.8:  (a)  Rate  of  heat  flow  at  the  inlet  and  exit  of  the  packed  bed 

(b)  Thermal  charging  of  the  packed  bed  (constant  temperature  walls} 
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3A.2.2  Effect  of  Darcv  number  (Da)  on  condensation 


The  effect  of  Darcy  number  on  condensation  was  also  studied  by  considering  three 
different  Da  values  for  a  fixed  Rep,  namely  1,500.  Different  Da  values  were  obtained  by 
changing  the  particle  size.  Figure  3.10  shows  the  variations  of  the  average  overall 
condensation  rate  and  the  total  condensate  accumulation  in  the  packed  bed  per  unit  width  as 
a  function  of  time  for  Da  values  of  1.49xl0'8,  2.66xl0'8  and  1.067xl0'7.  The 
corresponding  dp  values  were  0.75,  1.0  and  2.0  mm,  respectively.  As  dp  increases  and, 
hence,  Da  increases,  the  specific  surface  area  of  the  particles  decreases.  This  causes 
smaller  amount  of  heat  transfer  between  the  solid  and  fluid  phases  and  a  faster  advancement 
of  the  thermal  penetration  depth  of  the  fluid  phase  (faster  attainment  of  lower  densities  than 
the  saturation  vapor  density  at  the  initial  condensation  points).  Therefore,  although 
condensation  zone  advances  faster  the  condensation  rate  and  the  amount  of  condensate 
accumulation  remain  small.  Decreasing  the  dp  at  constant  Rep  (by  increasing  the  pressure 
difference  applied  across  the  packed  bed)  increases  the  vapor  velocities  and,  hence,  causes 
larger  mass  flow  rates  and  larger  condensation  rates.  Also  as  dp  is  reduced  the  specific 
surface  area  of  the  particles  increases.  The  heat  transfer  between  the  solid  and  fluid  phases 
becomes  more  vigorous,  energy  is  transferred  at  a  faster  rate  from  the  fluid  to  the  solid 
particles  and,  hence,  it  takes  longer  for  the  fluid  phase  at  the  initial  condensation  points  to 
reach  temperatures  high  enough  at  which  the  vapor  density  becomes  less  than  the  saturation 
vapor  density  and  condensation  stops.  Hence,  higher  condensation  rates  are  sustained  for 
longer  periods  of  time  at  these  points,  resulting  in  larger  condensate  accumulation. 
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Figure  3.10:  Effect  of  Darcy  number  on  condensation 


3.4.23  Effect  of  thermal  capacity  of  the  solid  phase  on  condensation 


In  older  to  analyze  whether  the  thermal  capacity  of  the  solid  phase  utilized  had  any 
effect  on  condensation,  a  diffeient  material,  namely  Iff  Carbon-steel,  was  employed  as  the 
solid  phase  material  for  qualitative  comparisons.  The  variations  of  the  average  overall 
condensation  tate  and  the  total  condensate  accumulation  per  unit  width  of  the  packed  bed  as 
functions  of  lime  are  depicted  in  Figure  3.11,  for  the  case  of  lead  and  steel  as  solid  material 
for  three  different  Rc*p  values.  As  may  be  seen  from  this  figure,  although  for  a  fixed  Re() 
the  condensation  rate  is  almost  the  same  at  the  beginning  for  both  solid  materials,  high 
condensation  rates  are  sustained  for  a  longer  period  of  time  in  the  case  of  steel  than  in  the 
case  of  lead.  The  icason  for  this  is  that  the  thermal  capacity  of  steel  per  unit  volume  is 
approximately  2.5  times  that  of  lead.  Due  to  this  fact  along  with  the  high  heat  transfer  rate 
between  the  solid  and  fluid  phases,  the  temperature  propagation  in  the  fluid  phase  will  be 
slower  in  the  case  of  steel  than  in  the  case  of  lead  because  more  of  the  thermal  energy  of  the 
fluid  is  transferred  to  the  solid  in  the  case  of  steel  than  in  the  case  of  lead.  In  the  case  of 
steel  the  vapor  phase  will  take  a  longer  time  to  reach  a  high  enough  temperature  at  which 
the  vapor  density  will  become  lower  than  the  saturation  vapor  density  corresponding  to  that 
temperature  and,  therefore,  the  condensation  durations  will  be  longer. 


3.5  Conclusions 

In  this  section,  the  development  of  a  model  for  analyzing  condensation  in  a 
convective  flow  ot  a  vapoi  through  a  packed  bed  was  carried  out.  This  was  achieved  by 
incorporating  the  basic  thetmodvnamic  fundamentals  of  phase  change  into  the  model. 


time  (sec) 
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Figure  3.11:  Effect  of  the  thermal  capacity  of  the  solid  phase  of  the  packed  bed  on 
condensation 
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From  the  results  of  the  case  studies  performed,  it  can  be  concluded  that  the  amount  of 
condensation  of  the  working  fluid  is  dependent  on  several  parameters.  For  a  Fixed  particle 
size,  the  higher  the  pressure  difference  applied  across  the  packed  bed,  i.e.,  higher  the 
nominal  Rep,  larger  will  be  the  rate  and  total  amount  of  condensation,  and  slightly  shorter 
will  be  the  duration  of  condensation  in  the  packed  bed.  This  is  a  consequence  of  the  fact 
that  the  working  fluid  has  larger  mass  flow  rates  and  is  under  larger  compression  forces. 

The  amount  of  condensation  in  the  packed  bed  was  also  found  to  be  influenced  by  the 
thermophysical  properties  of  the  bed  particles.  This  could  be  explained  by  the  heat 
interactions  between  the  working  fluid  and  the  solid  particles  which  determine  the  time 
taken  for  the  thermal  penetration  to  travel  in  both  the  solid  and  the  fluid  phases  in  the 
packed  bed.  Larger  thermal  capacitance  materials  for  the  bed  particles  will  result  in 
maintaining  the  temperature  of  the  vapor  low  and  the  density  of  the  vapor  high,  resulting  in 
a  larger  condensation  rate. 

Constant  temperature  wall  conditions  yielded  larger  amounts  of  condensation 
compared  to  insulated  wall  conditions,  proving  the  fact  that  continuous  heat  removal  from 
the  packed  bed  will  result  in  larger  amounts  of  condensation. 

Case  studies  in  Section  11  also  established  that  in  the  case  of  insulated  boundary 
conditions  with  immobile  liquid  phase  in  the  packed  bed,  the  application  of  one¬ 
dimensional  formulation  is  extremely  accurate  for  the  ranges  of  Rep  and  Da  considered  in 
this  study.  On  the  other  hand,  for  constant  temperature  boundary  conditions  two- 
dimensional  formulation  is  absolutely  necessary. 

For  the  t>pe  of  pioblem  considered  in  this  section,  the  amount  of  condensation  in  the 
packed  bed  is  \er>  sensitive  to  the  variation  of  the  fluid  temperature  on  which  the  saturation 
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vapor  density  is  highly  dependent.  Therefore,  it  is  very  crucial  to  formulate  the  problem 
with  no  local  thermal  equilibrium  assumption  between  the  solid  and  the  fluid  phases  in 
order  not  to  introduce  significant  inaccuracies  in  tracking  the  condensation  process. 


SECTION  IV 


ANALYSIS  OF  A  LATENT  HEAT  STORAGE  PACKED  BED 
UNDERGOING  A  CONDENSING  FLOW  OF  A  VAPOR 

In  this  section,  the  problem  investigated  in  Section  III  will  be  extended  in  complexity 
by  introducing  a  phase  change  material  (PCM)  to  replace  the  sensible  heat  storage  material 
of  the  bed  panicles.  Therefore,  in  addition  to  phase  change  in  the  working  fluid,  this 
problem  involves  phase  change  in  the  bed  panicles.  The  model  used  in  Section  III  is 
modified  such  that  it  will  accommodate  the  phase  change  processes  in  the  packed  bed 
panicles  for  the  present  problem.  The  numerical  solution  algorithm  that  shows  some 
appreciable  differences  with  that  of  the  previous  problem  is  discussed.  This  is  followed  by 
a  discussion  of  the  results  of  a  number  of  case  studies  performed. 


4.1  Statement  of  the  problem 

The  problem  under  consideration  is  the  thermal  charging  of  a  packed  bed  of  regularly 
sized  spheres  of  encapsulated  phase  change  material.  This  form  of  packed  bed  is  also 
known  as  a  latent  heat  storage  packed  bed.  The  schematic  diagram  of  the  problem 
considered  in  the  present  study  is  shown  in  Figure  4.1.  As  in  the  previous  section,  the 
packed  bed  is  initially  filled  with  R-12  vapor  at  a  slightly  superheated  state  at  uniform 
temperature  and  pressure.  Then  it  is  subjected  to  flow  of  superheated  R-12  from  a 
reservoir  which  has  higher  temperature  and  pressure  than  those  initially  prevailing  in  the 
packed  bed.  The  physical  conditions  for  this  problem  are  very  similar  to  those  of  the 
previous  section.  In  this  problem,  however,  the  packed  bed  incorporates  a  PCM  whose 
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Figure  4.1:  Schematic  diagram  of  the  problem 
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melting  temperature  falls  within  the  operating  temperature  range.  Therefore,  during  the 
course  of  the  thermal  charging  process  of  the  packed  bed  there  will  be  phase  change 
(melting)  within  the  encapsulated  PCM  particles. 


4.2  Mathematical  model 

The  volume  averaging  technique  employed  previously  is  also  utilized  in  the  present 
study.  All  of  the  eight  assumptions  and  simplifications  used  in  Section  3.2  are  assumed  to 
be  valid  here  too.  In  addition,  the  following  assumption  is  also  employed  in  the  analysis  of 
this  problem:  the  temperature  within  each  individual  bed  particle  is  assumed  to  be  uniform 
at  any  instant,  and  the  encapsulating  material  is  assumed  to  have  the  same  thermophysical 
properties  as  the  PCM. 

The  governing  equations  for  the  present  problem  are  similar  to  those  of  the  previous 
problem  except  for  the  solid  phase  energy  equation.  Where  a  similar  energy  equation  can 
be  used  when  the  PCM  is  not  undergoing  phase  change,  this  equation  becomes  inapplicable 
when  phase  change  is  taking  place  in  the  PCM.  During  the  process  of  phase  change  in  the 
PCM,  the  PCM  temperature  is  assumed  to  remain  constant  at  a  value  equal  to  its  melting 
temperature.  The  governing  equations  for  this  problem  can  be  summarized  as  follows: 

Vapor  phase  continuity  equation: 

J^(ey<  py>Y)  +  V.(<py>Y<  V)  =  ~<m>  (4.1) 

Liquid  phase  continuity  equation: 
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Vapor  phase  equation  of  motion  : 


V  <  PY  >T  =  - 


Y_  <Py>Tf£y 


<?y>  *iy 


1/2  [<7Y>-<^Y>]|<^Y>j  “  Ky<^> 


which  takes  the  following  one-dimensional  form: 


d  n  y  <Py>  F£y  2  Py 
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Liquid  phase  equation  of  motion: 
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Fluid  phase  energy  equation : 
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Solid  phase  (PCM)  energy  equation: 

When  there  is  no  phase  change  in  the  PCM  we  have: 


'i  * 


(4.2) 


(4.3) 


(4.3a) 


(4.4) 


] 


(4.5) 
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eoP<»(cpo  ~  ^•[k0cff^<To >0] 

-  h^a„p[<  T»  ><’-<  Tr  ‘  S  a«r[<  T»  >0- <  Tf >f  ] 


(4.6a) 


and  when  the  PCM  is  undergoing  phase  change  we  have: 


<T „>0  =  T  „ 

a  melt 


(4.6b) 


Volume  constraint  relation: 


£<,  +  £^(0  +  ^(0  =  1 


(4.7) 


Equation  of  state  for  vapor  phase: 


<Py>7=  <pY>YRY<Tf  >f 


(4.8) 


Thermodynamic  relation  for  the  saturation  density  of  vapor: 


exp  [  A 


Py.s " 


A| 
Tr  J 


RYTf 


(4.9) 


where  A,  and  B  are  known  constants,  Tf  is  in  degrees  Kelvin  and  P  y,  s  is  in  kg/m3. 


The  modeling  of  the  other  parameters  including  the  effective  thermal  conductivities  of 
the  solid  (PCM)  and  the  fluid  phases,  k^ff  and  kfcff,  the  permeability  of  the  packed  bed,  K, 
the  permeability  for  the  vapor  phase,  Ky,  the  geometiic  factor,  F,  the  relative  permeability 
of  the  liquid  phase,  krp,  the  empirical  correlations  for  the  fluid-to  particle  heat  transfer 
coefficient,  haj,  and  the  specific  surface  area  of  the  packed  bed  for  the  vapor  phase  and  for 
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die  liquid  phase,  a^y  and  a^,  is  carried  out  exactly  in  the  same  way  as  .was  done  in  Section 
III  in  equations  (3.10)  through  (3.18). 

4.2.1  Initial  and  boundary  conditions 

Following  the  arguments  in  Section  3.2.1,  the  initial  conditions  for  the  present 
problem  are  t  mathematically  given  as: 

Tf(x,y,t=0)  =  Ts(x,y,t=0)  =  T0 

Py(x,y,t=0)  =  P0  (4.10) 

Uy(x,y,t=0)  =  0 

There  is  a  continuous  flow  of  high  temperature  vapor  into  the  packed  bed  from  a 
reservoir  while  the  pressure  at  the  exit  of  the  packed  bed  is  maintained  at  the  initir.l  bed 
pressure.  The  following  mathematical  forms  express  these  boundary  conditions: 

Tf(x=0,y,t)  =  Tin 

Py(x=0,y,t)  =  Pin  at  t>0+  (4.11) 

Py(x=L,y,t)  =  Pout  =  P0 

where: 

T0  =  300  K  P0  =  796kPa  Tin  =  350K  Pout-796kPa  and  Pin  =  811.2kPa 


Since  we  were  primarily  interested  in  the  energy  storage  characteristics  of  the  packed 
beds,  insulated  wall  boundary  conditions  which  were  the  most  appropriate  ones  have  been 
used  in  this  section.  These  are  given  as: 


9Tf 
k  feff'aT' 


y  =  0,  y  =  H 


=  0 
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y  =  0,  y  =  H 


(4.12) 


4.2.2  Physical  conditions  for  the  numerical  runs 

In  the  present  problem  two  different  phase  change  materials  (latent  heat  storage 
materials)  and  a  sensible  heat  storage  material  were  considered  for  qualitative  comparisons 
in  the  thermal  charging  behavior  of  the  packed  bed.  The  PCMs  selected  for  this  study  were 
myristic  acid  which  has  a  melting  temperature  equal  to  331  K,  and  lithium-nitrate-trihydrate 
which  has  a  melting  temperature  of  303  K.  We  will  refer  to  these  materials  as  PCM1  and 
PCM2,  respectively.  The  sensible  heat  storage  material  chosen  was  1%  Carbon  steel. 

The  physical  property  values  of  the  materials  which  were  used  in  the  numerical 
computations  are  as  follows: 


R-12: 

ky  =  0.0097W/m.K  (cp)  =  71 07/kg. K  py  =  12.6xl0-6kg/m.s  Ry  =  0.0687588J/kg.K 
Ahvap  =  1 1 1300J/kg  (cp)p  =  1 1 15J/kg.K  =  179.2xl0-6kg/m.s  pp  =  1 190.35kg/m3 
kp  «  0.0545  W/m.K 

Myristic  acid  (PCM1): 


cp  =  1590J/kg.K 

k  =  0.1  W/m.K 

p  =  860kg/m3 

for  solid  phase 

cp  =  2260J/kg.K 

k  =  0.1  W/m.K 

p  =  860kg/m3 

for  liquid  phase 

Lithium-nitrate-trihvdrate  (PCM2): 

cp  =  2090J/kg.K 

k  =  0.5  W/m.K 

p  =  1550kg/m3 

for  solid  as  well  as  liquid  phase 

1%  Carbon-steel: 


Cp  =  473J/kg.K  k  =  43  W/m.K  p  =  7800kg/m3 


4.3  Solution  procedure  ' 

The  numerical  scheme  used  for  the  previous  problem  is  used  here  with  some 
modifications  that  accommodate  the  phase  change  processes  in  the  PCM.  These  include 
some  flags  in  the  numerical  code  that  perform  the  checks  for  differentiating  which  one  of 
equations  (4.6a)  and  (4.6b)  will  be  used  for  computing  the  PCM  temperature. 

In  the  solution  of  this  problem,  depending  on  whether  phase  change  is  occurring  in 
either  the  working  fluid  or  the  PCM  at  a  given  location  and  instant  in  the  packed  bed,  the 
governing  equations  and  the  solution  format  take  different  forms.  It  is  assumed  that  the 
condensation  occurs  when  the  density  of  the  vapor  becomes  equal  to  the  saturation  vapor 
density. 

When  there  is  no  phase  change  in  either  the  working  fluid  or  the  bed  particles,  the 
field  variables  <  Py  >Y,  £p,  <uy>,  <Vp>,  <Tf>f,  <T0>°,  and  <Py>Y 

are  determined  from  equations  (4.1)  to  (4.8),  respectively.  At  the  same  time  Py,  s  is 
determined  from  equation  (4.9).  This  is  necessary  for  determining  when  to  switch  to  the 
special  solution  format  used  for  situations  in  which  condensation  takes  place  in  the  working 
fluid.  It  is  carried  out  in  order  to  determine  when  to  switch  from  noncondensing  mode  to 
condensing  mode  solution  format. 

The  phase  change  process  in  the  PCM  also  needs  special  consideration.  The 
following  physical  characteristics  of  the  PCM  are  built  into  the  numerical  code.  Once  the 
PCM  at  a  certain  location  and  time  reaches  its  melting  temperature  during  the  thermal 
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charging  of  the  packed  bed,  its  temperature  remains  constant  until  phase  change  (melting)  is 
complete  in  the  capsules  at  that  location.  During  this  period,  the  solid  phase  (or  PCM) 
energy  equation  should  not  be  used  for  determining  the  PCM  temperature.  However,  the 
amount  of  heat  that  is  transferred  to  the  PCM  while  it  is  undergoing  melting  is  integrated  in 
time  in  order  to  determine  when  the  phase  change  is  completed.  Once  the  phase  change  is 
completed  in  the  PCM  at  a  certain  location,  the  procedure  of  determining  the  PCM 
temperature  is  switched  back  to  solution  from  the  regular  PCM  energy  equation,  i.e., 
equation  (4.6a),  with  appropriate  liquid  PCM  properties  incorporated  into  the  numerical 
code. 

At  each  time  step,  after  the  PCM  temperature  reaches  the  melting  temperature,  the 
convective  heat  transfer  rate  from  the  working  fluid  to  the  PCM  and  conduction  heat 
transfer  rates  to  and  from  the  PCM  are  summed  up  for  each  node  and  taken  into  account. 
Furthermore,  during  a  given  time  step  both  the  net  energy  input  to  the  PCM  and  the  net 
total  energy  accumulated  in  the  PCM  up  to  the  end  of  that  time  step  are  computed  on  a  unit 
volume  basis  for  each  node  where  the  PCM  is  undergoing  phase  change.  This  process  is 
carried  out  until  the  net  total  energy  accumulated  in  the  PCM  becomes  equal  to  the  latent 
heat  of  fusion  of  the  PCM.  During  this  period  the  solid  temperature  is  taken  to  be  equal  to 
T’melt  as  given  by  equation  (4.6b).  After  this  period  <  Tc  >°  is  again  computed  from 
equation  (4.6a). 

Implementation  of  the  stability  of  the  numerical  scheme  and  the  accuracy  of  the 
solutions  were  carried  out  as  explained  in  the  previous  sections.  A  41  x  21  grid 
configuration  (which  corresponds  to  a  dimensionless  Dx  or  Dy  of  0.025)  was  found  to 
yield  accurate  enough  solutions. 
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4.4  Results  and  discussions 


It  was  again  found  convenient  to  nondimensionalize  some  of  the  field  variables  in 
presenting  the  numerical  solutions,  and  keep  others  in  dimensional  form  to  have  a  better 

y  y 

understanding  of  the  variations  of  the  important  parameters.  Variables  <PY>  ,  <  pY  >  , 
and  <  UY  >  are  nondimen sionali zed  with  respect  to  the  corresponding  reference  quantities 
P*,  p*  and  v*.  However,  the  pertinent  parameters  of  the  thermal  charging  of  the  packed 
bed  and  the  condensation  of  the  working  fluid  are  presented  in  dimensional  form  per  unit 
width  of  the  packed  bed  as  a  function  of  the  dimensional  time.  The  value  chosen  for  P* 
was  100  kPa,  while  p*  was  then  computed  from  the  equation  of  state  using  P*  and  the 
initial  temperature,  T0.  The  calculation  of  the  reference  velocity,  v*,  was  based  on  the 
vapor  phase  momentum  equation  by  incorporating  the  global  pressure  difference  applied 
across  the  packed  bed  and  the  vapor  density  as  calculated  from  the  equation  of  state  by 
using  the  mean  value  of  the  inlet  and  exit  pressures. 

Due  to  the  physical  conditions  considered  in  the  present  study,  the  maximum  local 
liquid  fraction  of  the  working  fluid  did  not  exceed  the  critical  value  above  which  the  liquid 
becomes  mobile,  i.e.,  the  liquid  phase  was  always  in  pendular  state.  This,  combined  with 
the  fact  that  only  insulated  boundary  conditions  were  considered  in  the  present  study,  made 
the  problem  essentially  one-dimensional.  It  was  checked  through  numerical 
experimentation  that  the  solution  of  the  one-dimensional  form  of  the  governing  equations 
did  not  have  any  appreciable  difference  from  the  solution  to  the  two-dimensional  form. 
Therefore,  the  results  of  the  one-dimensional  solution  will  be  presented  for  convenience. 

Two  distinct  stages  that  were  observed  in  the  previous  two  problems  can  be  easily 
identified  in  the  solution  of  the  present  problem  too.  The  early  stage  with  very  strong 
transient  effects,  i.e.,  drastic  spatial  and  temporal  changes  in  the  field  variables,  lasts  for  a 
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very  short  time.  During  the  early  stage  the  pressure  distribution  across  the  packed  bed 
evolves  and  assumes  an  almost  linear  form  which  is  maintained  afterwards  during  the  rest 
of  the  charging  process  which  will  again  be  referred  to  as  the  later  stage . 

The  first  PCM  employed  for  the  encapsulated  bed  particles  is  myristic  acid  (or 
PCM1),  which  has  a  melting  temperature  of  331  K.  The  results  for  the  case  in  which  this 
material  is  used  will  be  presented  in  detail.  Once  the  high  pressure,  high  temperature  vapor 
is  applied  at  the  inlet  of  the  pack  d  bed,  the  vapor  moving  through  the  packed  bed  becomes 
compressed.  Since  it  gives  most  of  its  excess  internal  energy  to  the  colder  PCM  particles, 
its  density  reaches  the  saturation  vapor  density  at  certain  locations  and  condensation  takes 
place.  Most  of  this  condensation  occurs  in  the  early  siage  while  the  pressure  distribution 
linearizes  and  the  density  of  the  vapor  adjusts  itself  accordingly.  Afterwards,  the  vapor 
reaches  superheated  conditions  at  all  locations  and  no  more  condensation  takes  place. 

The  early  stage  variations  of  the  density,  velocity  and  the  pressure  of  the  working 
fluid,  and  the  liquid  fraction  are  shown  in  Figure  4.2.  The  early  stage  is  somewhat 
extended  such  that  it  includes  the  period  during  which  more  than  99%  of  the  condensation 
in  the  working  fluid  is  completed.  The  consequence  of  the  step  change  boundary 
conditions  can  easily  be  seen  from  the  high  velocities  at  the  inlet  of  the  packed  bed.  These 
high  velocities  die  out  as  the  pressure  distribution  becomes  linear.  Except  for  a  short 
thermal  entry  region,  during  the  early  stage  the  changes  in  the  field  variables  are  mostly 
pressure  dependent  because  there  is  no  appreciable  thermal  penetration.  It  should  be  noted 
that  since  the  liquid  fraction  never  reaches  the  critical  value  for  becoming  mobile,  the  eb 

distribution  remains  the  same  throughout  the  later  stage.  The  early  stage  results  for  this 
problem  are  qualitatively  similar  to  the  ones  of  the  problem  in  Section  III. 
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Figure  4.2:  Variations  of  different  field  variable  distributions  during  the  early  stage 


The  changes  in  the  field  variables  during  the  later,  stage  can  be  attributed  to  die 
development  of  the  thermal  penetration  depth.  As  may  be  seen  from  Figure  4.3,  the 
pressure  distribution  remains  linear.  The  variations  of  the  temperatures  of  the  working 
fluid  and  the  PCM  are  very  smooth  until  the  PCM  reaches  its  melting  temperature.  After 
the  onset  of  melting  in  the  PCM,  a  distinct  discontinuity  can  be  observed  in  the  smoothness 
of  the  PCM  temperature  distribution.  This  is  because  for  a  certain  length  of  the  packed  bed 
there  is  no  change  in  the  PCM  temperature  while  the  PCM  is  undergoing  melting.  In  this 
region  the  working  fluid  also  adjusts  itself  accordingly.  This  can  be  observed  in  the 
temperature  profiles  for  time  levels  2  through  8  in  Figure  4.3(a).  This  kind  of  qualitative 
behavior  can  also  be  seen  in  the  vapor  velocity  and  density  variations  along  the  packed  bed 
at  the  corresponding  time  levels.  This  is  a  consequence  of  the  fact  that  the  vapor  density 
variation  is  precisely  related  to  the  pressure  and  temperature  variations.  It  can  also  be 
noticed  that  there  is  an  inverse  relationship  between  the  vapor  density  and  vapor  velocity. 
This  is  due  to  the  fact  that  the  transient  term  in  the  vapor  continuity  equation  loses  its 
dominance  and  the  convective  term  dominates  during  the  later  stage.  Once  the  packed  bed 
becomes  thermally  charged,  then  the  vapor  density  distribution  becomes  linear  similar  to 
the  pressure  distribution  as  can  be  explained  by  the  equation  of  state. 

The  variation  of  the  average  overall  condensation  rate  and  the  accumulative 
condensate  per  unit  width  of  the  packed  bed  are  shown  in  Figure  4.4.  The  overall 
condensation  rate  was  computed  by  integrating  the  condensation  rates  at  the  individual  grid 
points  over  the  associated  volumes.  Integration  of  this  over  time  yielded  the  accumulative 
condensate  per  unit  width  of  the  packed  bed.  Again  the  qualitative  behavior  in  these  figures 
is  very  similar  to  the  ones  in  the  previous  section. 
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Of  interest  to  the  thermal  charging  of  the  packed  bed  are  the  rates  of  heat  flowing 
into  and  out  of  the  packed  bed  as  a  function  of  time.  These  were  also  computed  for  the 
given  crosssection  of  the  packed  bed  with  a  unit  width.  The  variations  of  these- 
quantities  are  depicted  in  Figure  4.5(a)  except  for  a  very  short  time  section  at  the  beginning 
of  the  charging  process,  which  was  left  out  for  obtaining  a  better  scale  on  the  figure.  The 
variation  of  the  net  energy  stored  per  unit  width  of  the  packed  bed  as  a  function  of  time  is 
also  depicted  in  Figure  4.5(b).  At  the  beginning  of  the  charging  process,  the  vapor  flowing 
out  of  the  packed  bed  leaves  at  a  low  temperature  close  to  the  initial  temperature  and,  thus, 
there  is  a  large  difference  between  the  heat  flow  rates  at  the  inlet  and  at  the 'exit  of  the 
packed  bed.  Hence,  the  rate  of  energy  storage  is  large.  Afterwards,  for  a  major  portion  of 
the  charging  process  there  is  a  uniform  difference  between  the  rates  of  heat  flowing  into 
and  out  of  the  packed  bed  causing  a  linear  increase  in  the  amount  of  energy  stored.  Once 
the  phase  change  is  complete  in  all  particles  of  the  packed  bed,  both  the  working  fluid 
temperature  and  the  PCM  temperature  at  the  exit  of  the  bed  rise  rapidly  causing  a  rapid 
decrease  in  the  gap  between  the  heat  flow  rates  at  the  inlet  and  the  exit.  When  the  packed 
bed  becomes  completely  charged  thermally,  there  remains  no  difference  between  the 
amount  of  heat  flowing  into  and  out  of  the  packed  bed. 


4.4.1 


Besides  myristic  acid  (PCM  1 ),  two  more  materials  were  considered  for  the  particles 
of  the  packed  bed.  These  were  lithium-nitrate-trihydrate,  or  PCM2,  and  1%  Carbon-steel. 
Computational  runs  were  made  for  these  cases  with  the  same  initial  and  boundary 
conditions  as  in  the  case  of  PCM1.  Figure  4.6  depicts  the  overall  average  condensation 
rate  and  condensate  accumulation  histories  for  all  three  cases.  Since  condensation  in  the 


Accumulative  condensate  (kg)  Overall  condensation  rate  (kg/sec) 


working  fluid  takes  place  in  a  very  short  span  of  time  at  the  beginning,  with  no  significant 
thermal  penetration  in  the  packed  bed,  the  ditference  in  the  results  for  the  three  materials 
considered  can  be  attributed  mainly  to  their  physical  properties.  An  analysis  shows  that 
the  thermal  capacitance  per  unit  volume  of  PCM2  is  only  approximately  12%  different  from 
that  of  steel,  where  the  thermal  capacitance  per  unit  volume  of  PCM1  is  almost  3  times 
smaller  than  that  of  steel.  Due  to  this  fact  and  the  high  heat  transfer  rate  from  the  working 
fluid  to  the  solid  phase,  the  temperature  propagation  in  both  the  solid  phase  as  well  as  the 
working  fluid  is  slower  in  the  cases  with  steel  and  PCM2  than  in  the  case  with  PCM1. 
Thus,  for  the  former  two  cases  it  takes  longer  time  for  condensation  to  end.  This  is  as  a 
result  of  the  longer  time  needed  for  the  vapor  to  reach  at  a  high  enough  temperature  at 
which  the  vapor  density  becomes  lower  than  the  saturation  vapor  density.  Consequently, 
higher  condensation  rates  are  sustained  for  longer  periods  and  larger  condensate 
accumulation  takes  place  in  the  case  of  steel  and  PCM2. 

4.4.2  Qualitative  comparison  of  the  thermal  charging  process 

The  three  materials  considered  for  [red  particles  in  this  section  are  the  same  as  those  of 
the  previous  section.  Time  histories  of  the  rates  of  heat  flow  into  and  out  of  the  packed  bed 
for  the  three  cases  are  shown  in  Figure  4.7,  while  that  of  the  net  energy  stored  in  the 
packed  bed  per  unit  width  is  shown  in  Figure  4.8.  It  can  be  seen  that,  although  at  the 
beginning  of  the  charging  process  the  energy  storage  is  almost  the  same  for  all  three  cases, 
it  shows  a  different  variation  during  later  times.  While  that  of  the  case  with  steel  looks  like 
a  conventional  charging  curve  which  decays  exponentially  with  time,  the  cases  with  the 
PCMs  have  a  linear  variation  for  a  major  portion  of  the  charging  period.  These  linear 
portions  correspond  to  the  time  spans  during  which  the  temperature  of  the  working  fluid  at 
the  exit  of  the  packed  bed  is  fairly  constant  and  approximately  equal  to  the  melting 
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Figure  4.7:  Heat  flow  rate  into  and  out  of  the  packed  bed  for  cases  with  different  bed  particle  materials 


temperature  of  the  PCM,  since  once  the  PCM  reaches  it  melting  temperature  its  temperature 
remains  constant  until  phase  change  is  complete  in  the  PCM.  During  this  period  the  vapor 
temperature  cannot  drop  below  the  melting  temperature  of  the  PCM  and,  thus,  we  have  a 
fairly  constant  vapor  exit  temperature.  In  order  to  illustrate  this  more  clearly,  the  time 
history  of  the  vapor  exit  temperature  is  presented  in  Figure  4.9  for  the  three  cases.  As  may 
be  seen  from  this  figure,  the  vapor  exits  the  packed  bed  at  a  lower  temperature  for  a  longer 
period  of  time  in  the  case  of  PCM2  than  in  the  case  of  PCM1.  This  is  partly  due  to  the 
lower  melting  temperature  of  PCM2  and  partly  due  to  its  larger  thermal  capacitance. 
Therefore,  the  difference  between  the  rates  of  heat  flowing  into  and  out  of  the  bed  is  larger 
and  lasts  longer  in  this  case,  thus,  giving  a  larger  amount  of  net  energy  stored. 

Figure  4.8  also  shows  that  using  a  PCM  with  a  certain  melting  temperature  may  not 
always  be  a  better  choice  over  a  sensible  heat  storage  material.  For  instance,  for  the 
boundary  conditions  and  the  size  considered  in  the  present  Work,  steel  seems  to  perform 
better  than  PCM1  if  we  are  interested  in  the  energy  storage  range  of  up  to  2200  kJ  per  unit 
depth  over  approximately  50  seconds  with  no  limitation  on  the  weight  of  the  packed  bed. 


4.5  Conclusions 

In  this  section,  in  which  the  emphasis  was  given  to  the  study  of  the  behavior  of 
energy  storage  packed  beds,  it  was  shown  that  the  characteristics  of  the  thermal  charging  of 
the  latent  heat  storage  packed  beds  differed  appreciably  from  those  of  the  sensible  heat 
storage  packed  beds.  While  the  thermal  charging  curve  of  a  sensible  heat  storage  packed 
bed  looks  very  much  like  the  conventional  charging  curve  of  a  capacitor,  that  of  a  latent 
heat  storage  packed  bed  shows  different  behavior  in  that  its  variation  is  linear  in  time  during 
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Figure  4.9:  Time  history  of  vapor  exit  temperature  for  cases  with  different  bed  particle  materials 


a  major  part  of  the  charging  process.  This  behavior,  which  is  the  consequence  of  a  fairly 
constant  vapor  exit  temperature  during  the  part  of  the  charging  process  in  which  phase 
change  takes  place  in  the  PCM,  was  found  to  be  dependent  on  the  melting  temperature  of 
the  PCM.  Case  studies  performed  showed  that  energy  storage  in  a  packed  bed  employing  a 
PCM  would  be  more  efficient  if  the  melting  temperature  of  the  PCM  is  close  to  the  lower 
end  of  the  operating  range  of  the  packed  bed.  Because  this  way  optimum  use  of  the  storage 
capacity  would  be  realized,  i.e„  the  least  amount  of  energy  would  be  able  to  escape  from 
being  stored  in  the  bed  particles. 

In  the  investigations  in  Section  IV  the  differences  between  the  solid  and  fluid 
temperatures  became  even  more  apparent  than  in  the  previous  cases  in  Sections  II  and  III, 
confirming  the  necessity  of  modeling  such  problems  with  no  local  thermal  equilibrium 
assumption. 
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SECTION  V 


ANALYSIS  OF  THERMAL  CHARGING  AND  DISCHARGING 
OF  SENSIBLE  HEAT  AND  LATENT  HEAT 
STORAGE  PACKED  BEDS 


5.1  Introduction 

The  operation  of  the  packed  bed  which  employs  encapsulated  PCM  as  energy 
storage  medium  is  of  interest  for  the  operation  during  thermal  chaiging  and  discharging 
modes.  It  is  important  to  analyze  the  transient  response  of  the  packed  bed  during  these 
modes  of  operation  in  order  to  deteimine  the  duration  for  energy  storage  or  removal,  or  for 
determining  the  size  of  the  packed  bed  required  for  a  given  operating  range.  In  this  section, 
thermal  charging  and  discharging  of  sensible  heat  and  latent  heat  storage  packed  beds  are 
modeled  and  simulated.  The  results  of  this  investigation  show  that  the  energy  storage 
characteristics  of  a  latent  heat  storage  packed  bed  are  quite  different  than  those  of  a  sensible 
heat  storage  packed  bed,  and  that  they  are  very  much  dependent  on  the  melting/freezing 
temperature  of  the  phase  change  material  utilized  for  energy  storage. 


5.2  Problem  statement  and  formulation 

The  physical  system  considered  consists  of  a  horizontal  channel  filled  with 
randomly  packed  fixed  particles  of  regularly  sized  and  shaped  spheres,  forming  a  packed 
bed  as  shown  in  Figure  5. 1 .  The  depth  of  the  packed  bed  is  assumed  to  be  large  enough  in 
order  to  avoid  three-dimensionality  of  the  problem.  Initially  the  void  volume  of  the  packed 
bed  is  filled  with  quiescent  working  fluid  which  is  at  uniform  temperature  and  pressure  and 
in  thermal  equilibrium  with  the  bed  particles.  Vapor  from  a  reservoir  at  a  higher 
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Figure  5. 1 :  Schematic  diagram  of  the  problem. 
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temperature  and  pressure  is  allowed  to  flow  through  the  packed  bed  in  the  thermal  charging 
mode.  In  the  thermal  discharging  mode,  vapor  at  a  lower  temperature  is  considered  to  flow 
through  the  packed  bed.  Relatively  high  particle  Reynolds  number  flows  are  considered  in 
this  work,  assuring  that  the  flow  of  the  Refrigerant- 12  vapor  is  forced  convective  flow  in 
nature. 

t 

The  governing  equations  for  the  present  problem  were  developed  by  use  of  the  well 
known  "local  volume  averaging"  technique  and  they  follow  from  the  previous  work  of 
Vafai  and  Sozen  (1990).  It  should  be  emphasized  that  in  modeling  the  physical 
phenomena,  the  thermophysical  properties  of  the  encapsulation  material  were  assumed  to 
be  essentially  the  same  as  those  of  the  PCM  in  the  case  of  latent  heat  storage  packed  bed. 
Due  to  the  insulated  boundary  conditions  employed  in  this  investigation,  the  problem 
essentially  reduces  to  be  one-dimensional.  The  governing  equations  are  given  as  follows: 

Vapor  phase  continuity  equation: 

^(«<P.>')+£(<P.  >•<«»>)=  o  (5.1) 


Vapor  phase  momentum  equation: 
d<Pv>v  <  pv  >v  Fe 


dx 


K 


1/2 


<  «v  >2  "  YK  U *  > 


(5.2) 


Vapor  phase  energy  equation: 


£<  Pv>v  CP 


d<Tv>v 

dt 


+  Cp<pv  >V<MV>- 


d<T>y 


dx 


=  JxWn  +  K  {<  Ts  >'  -  <  Tv  >'} 


(5.3) 
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Solid  phase  (or  PCM)  energy  equation  : 


„  x  d<Ts>s  d[,  d<Ts>* 


-  <  t,  >'} 


Vapor  phase  equation  of  state: 

cPv>v  =<pv>vR<Tv>v  (5.5) 

where  subscript  or  superscript  "s"  denotes  the  solid  phase  in  the  case  of  sensible  heat 
storage  packed  bed  or  the  PCM  in  the  case  of  latent  heat  storage  packed  bed. 

The  above  set  of  equations  forms  a  system  of  five  equations  in  five  unknowns, 
namely,  <  pv  >v,  <  uv  >,  <  Tv  >v,  <  Ts  >s,  and  <  Pv  >v. 

The  initial  conditions  of  the  problem  were  as  follows: 

Pv  (x,  t  =  0)  =  P0 

Tv(x,t  =  0)  =  Ts(x,t  =  0)  =  T0 

uv(x,t  =  0)  =  0  (5-6) 

and  the  corresponding  initial  values  of  py  are  calculated  from  eqn.  (5). 

The  boundary  conditions  of  the  problem  were  as  follows: 

Pv(x  =  0,t)  =  Pm 
Pv(x  =  L,t)  =  P0 

T,(x  =  0,t)  =  T,„  (3J) 

where  P0  =  lOOkPa,  Pm  =  106.83kPa,  T0  =  300  K  and  TVU1  =  350  K  for  thermal  charging 
mode,  and  T™  =  300  K  for  thermal  discharging  mode.  The  nominal  particle  Reynolds 
number,  Rep)  for  these  boundary  conditions  was  taken  as  1,000. 
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Modeling  of  the  permeability  of  the  packed  bed,  Kv,  and  geometric  function  Em 
the  vapor  phase  momentum  equation,  the  specific  interfacial  surface  area  of  the  packed  bed, 
the  fluid-to-particle  heat  transfer  coefficient  and  the  effective  thermal  conductivities  of  the 
two  phases  were  carried  out  exactly  in  the  same  way  as  in  Section  II.  Also,  an  average 
value  of  0.39  was  used  for  the  porosity  of  the  packed  bed  based  on  the  findings  of 
Benanati  and  Brosilow  (1962).  The  diameter  of  the  bed  particles  was  taken  to  be  2  mm. 


The  thermophysical  properties  used  in  the  present  investigation  were  as  follows: 


Refrigerant- 12 

cp  =  7 10  J  kg'1  K"1 
k  =  0.0097Wm“'K"1 
p.  =  12.6x10**  kgm"1  s"1 
R  =  68.7588  J  kg'1  K"1 


1  %  Carbon-steel 


Mvristic  acid  (PCM) 


cp  =  473  J  kg'1  K"1  cp  =  1590Jkg-1  K"1  (solid) 

k  =  43  Wm'1  K"1  cp  =  2260  J  kg"1  K"1  (liquid) 
p  =  7800  kg  m"3  k  =  0.1  Wm-1  K"1 
p  =  860kgm~3 
hIf  =  200.5xl03  Jkg"1 


5.3  Solution.  Emmlms 

It  should  be  noted  that  for  the  solution  of  the  problem  for  latent  heat  storage  packed 
bed  (LHSPB),  the  solution  algorithm  is  the  same  as  in  the  case  of  sensible  heat  storage 
packed  bed  presented  in  Section  2.3  except  that  during  the  thermal  charging  process,  once 
the  temperature  of  the  PCM  reaches  its  melting  temperature  it  remains  constant  until  melting 
of  the  PCM  is  complete.  This  condition  has  been  incorporated  into  the  numerical  algorithm 
such  that  at  a  given  point  once  the  melting  temperature  is  reached,  the  PCM  temperature  is 
kept  constant  while  the  net  rate  of  heat  flowing  into  the  associated  volume  is  continuously 
monitored  and  integrated  over  time  to  give  the  energy  that  is  used  to  melt  the  PCM  until  it 
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becomes  equal  to  the  amount  required  for  melting.  After  this  stage  the  PCM  temperature  is 
again  computed  from  eqn.  (5.4),  this  time  with  the  liquid  properties  of  the  PCM 
incorporated  in  eqn.  (5.4).  Similar  provisions  are  made  for  the  thermal  discharging 
process  of  the  LHSPB  which  involves  freezing  of  the  PCM. 


5-4  Results,  an, cl  .discussions 

The  two  problems  considered  in  the  present  work  were  the  complete  thermal 
charging  and  thermal  discharging  of  two  types  of  packed  beds,  one  with  sensible  heat 
storage  material,  chosen  to  be  1  %  Carbon-steel,  and  the  other  one  with  latent  heat  storage 
material,  chosen  to  be  myristic  acid  (a  PCM  with  a  melting  point  of  331  K).  Thermal 
charging  or  discharging  was  assumed  to  be  complete  when  the  difference  between  the 
incoming  and  exiting  vapor  temperatures  were  less  than  0.01  K. 

In  the  solution  of  the  thermal  charging  problem  it  was  found  that  two  distinct  stages 
were  present;  namely,  the  early  stage  and  the  later  stage.  These  have  been  discussed  in 
detail  in  the  previous  sections.  In  this  work,  emphasis  is  given  to  the  thermal  charging  and 
discharging  characteristics  of  the  two  types  of  packed  beds  using  a  compressible  working 
fluid. 


Figure  5.2  depicts  the  variation  in  the  temperature  distribution  of  the  vapor  and 
solid  phases  during  the  thermal  charging  processes  for  the  case  of  sensible  heat  storage 
material,  while  Figure  5.3  shows  the  corresponding  variations  in  the  case  of  the  latent  heat 
storage  material  (PCM).  The  results  presented  in  Figures  5.2  and  5.3  were  obtained  under 
identical  operating  conditions  (which  were  given  after  eqn.  (5.7)).  From  Figure  5.3  it  can 
be  easily  seen  how  the  thermal  charging  of  the  packed  bed  starts  with  sensible  heating  of 
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Figure  5.3:  Temperature  distribution  in  the  latent  heat  storage  packed  bed  during  thermal 
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the  PCM  particles  (as  can  be  seen  from  the  first  three  time  levels  presented)  and  once  the 
melting  temperature,  331  K,  of  the  PCM  is  attained,  the  PCM  temperature  remains  constant 
until  phase  change  is  complete  at  that  location  and  afterwards  sensible  heating  takes  place 
again.  Figures  5.4  and  5.5  depict  the  variation  of  the  corresponding  temperature 
distributions  for  the  thermal  discharging  modes  of  the  two  packed  beds  considered.  From 
these  figures  it  can  be  seen  that  in  the  case  of  sensible  heating  or  cooling  of  the  packed  bed, 
the  temperature  distribution  is  quite  smooth.  This  is  seen  both  in  Figures  5.2  and  5.4  as 
well  as  in  certain  portions  of  Figures  5.3  and  5.5.  In  Figure  5.3,  sensible  heating  of  a 
major  portion  of  the  packed  bed  occurs  during  the  earlier  stages  of  the  thermal  charging  as 
discussed.  Similarly,  it  can  be  seen  from  Figure  5.5  that  the  major  portion  of  the  sensible 
cooling  is  again  at  the  early  stages  of  the  thermal  discharging  (time  levels  1  through  5), 
during  which  the  temperature  drops  to  the  freezing  point  of  the  PCM. 

The  rate  of  heat  flow  (enthalpy  flux)  into  and  out  of  the  packed  bed  per  unit  width 
of  the  packed  bed  is  shown  in  Figure  5.6  for  the  thermal  charging  of  the  sensible  heat 
storage  packed  bed  (SHSPB).  The  decrease  in  the  heat  flux  into  the  SHSPB  is  due  to  a 
decrease  in  the  mass  flow  rate  of  the  working  fluid,  R-12.  This  was  an  interesting  finding 
for  the  compressible  flow  with  inertia  effects,  and  this  behavior  differed  completely  from 
an  incompressible  flow  modeled  by  Darcy's  law.  The  corresponding  decrease  in  the  mass 
flow  rate  at  the  exit  of  the  SHSPB  also  results  in  a  decrease  in  the  rate  of  heat  flow  out  of 
the  bed  at  the  beginning.  However,  as  the  temperature  of  the  bed  particles  and,  hence,  that 
of  the  working  fluid  at  the  exit  of  the  bed  start  increasing,  i.e.,  less  of  the  energy  of  the 
working  fluid  is  being  dumped  into  the  bed  particles,  the  heat  flux  out  of  the  SHSPB  starts 
increasing  until  the  inlet  and  exit  vapor  temperatures  become  equal  and  the  heat  fluxes  into 
and  out  of  the  bed  level  off.  It  should  be  noted  that  the  heat  flux  at  the  inlet  of  the  SHSPB 
starts  at  a  value  of  zero  at  time  equal  to  zero  and  sharply  increases  as  the  flow  of  hot  vapor 
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Figure  5.4:  Temperature  distribution  in  the  sensible  heat  storage  packed  bed  during  thermal 
discharging. 
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Figure  5.5:  Temperature  distribution  in  the  latent  heat  storage  packed  bed  during  thermal 
discharging. 
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info  the  bed  starts.  This  very  short  portion  of  the  charging  history  was  otnitted  from  Figure 
5.6  in  order  to  obtain  a  better  scald  for  the  Figure,  figure  5.7  shows  the  variation  of  heat 
fluxes  into  and  out  of  the  SHSPB  during  the  thermal  discharging  mode.  As  may  be  seen, 
due  to  the  decrease  in  the  vapor  inlet  temperature  (300  K)  there  is  a  sharp  decrease  in  the 
heat  flux  into  the  bed  at  the  beginning.  Later  increase  in  the  heat  flux  into  the  bed  is  due  to 
an  increase  in  the  mass  flow  rate  of  the  vapor  contrasting  to  the  decrease  during  the 
charging  process.  Earlier  increase  in  the  heat  flux  out  of  the  bed  is  due  to  the  same  reason, 
where  later  decrease  and  leveling  off  with  the  heat  flux  into  the  bed  can  be  explained  due  to 
the  fact  that  the  exit  temperature  of  the  vapor  phase  starts  decreasing  and  gradually  less 
amount  of  thermal  energy  is  being  removed  from  the  SHSPB.  The  time  history  of  the 
thermal  charging  and  discharging  of  the  SHSPB  is  shown  in  Figure  5.8.  This  Figure 
depicts  the  variation  of  the  net  energy  stored  per  unit  width  of  the  SHSPB.  In  this  case  the 
thermal  charging  and  discharging  times  of  the  bed  are  practically  equal. 

Figure  5.9  depicts  the  variation  of  the  heat  flux  in  continuous  form  for  the  complete 
thermal  charging  followed  by  the  complete  discharging  of  the  SHSPB.  The  corresponding 
variations  for  the  latent  heat  storage  packed  bed  (LHSPB)  are  depicted  in  Figure  5.10.  A 
few  interesting  points  can  be  observed  in  the  thermal  charging  and  discharging  of  the 
LHSPB.  Both  during  the  charging  and  the  discharging  modes  there  is  a  major  portion  of 
each  process  during  which  there  is  quite  a  uniform  difference  between  the  heat  fluxes  into 
and  out  of  the  LHSPB  which  is  not  seen  in  the  case  of  the  SHSPB.  The  reason  for  this 
kind  of  behavior  is  the  fact  that  during  these  portions  of  the  charging  or  the  discharging,  the 
temperature  of  the  PCM  in  the  particles  downstream  of  the  packed  bed  remains  constant  at 
the  melting/freezing  temperature.  Hence,  for  these  portions  of  the  charging  or  discharging 
processes,  the  vapor  exit  temperature  is  practically  equal  to  the  melting/freezing  temperature 
of  the  PCM  and,  therefore,  the  difference  between  the  inlet  and  exit  temperatures  for  the 
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Heat  flow  rate  into  and  out  of  the  SHSPB  during  thermal  discharging. 
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Figure  5.9:  Heat  flow  rate  into  and  out  of  the  SHSPB  during  complete  thermal  charging 

followed  by  complete  discharging. 
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followed  by  complete  discharging. 


vapor  phase  remains  practically  unchanged.  This  can  be  easily  seen  from  Figure  5.11 
which  depicts  the  time  history  of  the  vapor  exit  temperature  for  the  complete  thermal 
charging  and  discharging  of  both  the  SHSPB  and  the  LHSPB.  Another  interesting 
observation  can  be  made  with  respect  to  the  charging  and  discharging  times  of  the  LHSPB. 
From  both  Figures  5.10  and  5.1 1  it  can  be  easily  observed  that  the  thermal  discharging  of 
the  LHSPB  lasts  considerably  shorter  than  the  thermal  charging.  The  reason  for  this  can  be 
seen  from  these  figures  again.  For  the  major  portion  of  the  discharging  process  during 
which  the  temperature  of  the  PCM  downstream  of  the  LHSPB  remains  at  the 
melting/freezing  point,  there  is  a  constant  temperature  difference  between  the  vapor  inlet 
and  and  vapor  exit  temperatures.  This  difference  is  larger  than  the  corresponding 
difference  in  the  charging  process,  i.e.,  the  difference  between  the  vapor  inlet  temperature 
and  the  melting  point  of  the  PCM  during  the  charging  process  is  only  19  K  where  during 
the  discharging  process  the  difference  between  the  vapor  inlet  temperature  and  the  freezing 
point  of  the  PCM  is  3 1  K,  and,  hence,  energy  is  stored  at  a  lower  rate  in  the  charging 
process  than  energy  is  removed  during  the  discharging  process  for  the  portions  of  these 
processes  under  discussion. 

Figure  5.12  depicts  the  time  histories  of  the  net  energy  stored  in  the  packed  bed 
during  the  charging  and  the  discharging  modes  of  the  SHSPB  and  the  LHSPB 
respectively.  While  the  rate  of  energy  storage  and  removal  follows  a  similar  trend  in  the 
case  of  the  SHSPB,  i.e.,  very  fast  at  the  beginning  and  decaying  afterwards,  this  was  not 
the  behavior  in  the  case  of  the  LHSPB.  Again,  this  was  due  to  the  fact  that  a  major  portion 
of  the  packed  bed  remains  at  constant  temperature  (melting/freezing  point)  during  the 
thermal  charging/discharging  of  the  LHSPB. 
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Figure  5.11: 
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Figure  5. 12:  Thermal  charging  and  discharging  of  the  SHSPB  and  the  LHSPB. 


5.5  Conclusions 


The  transient  processes  of  thermal  charging  and  discharging  of  an  SHSPB  and  an 
LHSPB  with  a  compressible  working  fluid  were  simulated.  The  investigations  showed 
distinctly  different  energy  storage  characteristics  for  these  two  kinds  of  packed  beds.  The 
high  energy  storage  density  of  the  LHSPB  was  clearly  observed  from  the  studies  carried 
out.  For  the  two  energy  storage  materials  considered,  although  the  density  of  the  sensible 
heat  storage  material  (1%  Carbon  steel)  was  approximately  9  times  larger  than  that  of  the 
PCM  (myristic  acid),  the  total  energy  storage  capacity  of  the  LHSPB  was  higher. 

Use  of  a  compressible  working  fluid  and  accounting  for  the  inertia  effects  in  the 
vapor  phase  momentum  equation  resulted  in  a  time-dependent  mass-flow  rate  through  the 
packed  bed.  Utilization  of  non-local-thermal-equilibrium  analysis  was  found  to  be  crucial 
especially  in  the  case  of  LHSPB. 

Also,  in  the  case  of  the  LHSPB  it  was  observed  that  the  thermal  charging  and 
discharging  times  differed  considerably.  The  main  reason  for  this  difference  is  that  for  a 
major  portion  of  the  charging  or  discharging  process  the  PCM  temperature  downstream  of 
the  packed  bed  remains  constant  at  the  melting/freezing  temperature  of  the  PCM.  An 
important  conclusion  which  was  obtained  from  our  results  is  that  closer  the  PCM 
melting/freezing  temperature  to  the  charging  temperature  (vapor  inlet  temperature  during 
thermal  charging  mode)  longer  will  be  the  time  taken  for  charging  and  shorter  will  be  the 
time  taken  for  discharging  of  the  packed  bed.  Likewise,  closer  the  PCM  melting/freezing 
temperature  to  the  discharging  temperature  (vapor  inlet  temperature  during  thermal 
discharging  mode)  longer  will  be  the  time  taken  for  discharging  and  shorter  will  be  the  time 
taken  for  charging  of  the  packed  bed. 
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SECTION  VI 

ANALYSIS  OF  OSCILLATING  COMPRESSIBLE  FLOW 
THROUGH  A  PACKED  BED 


6.1  Introduction 

In  storage  of  thermal  energy  in  packed  beds,  the  dynamic  behavior  of  the  packed 
bed  system  is  an  important  consideration.  Ideal  conditions  of  constant  temperature  and 
pressure  inlet  conditions  are  quite  difficult  to  maintain.  Therefore,  examining  those 
conditions  more  closely  approximating  the  real-life  situations,  such  as  variable  pressure  or 
variable  temperature  inlet  conditions,  provides  better  insight  to  such  problems. 

Rigorous  models  have  been  developed  by  Vafai  and  Sozen  (1990),  and  Sozen  and 
Vafai  (1990)  for  the  forced  convective  flow  of  a  superheated  ideal  gas  and  forced 
convective  condensing  flow  of  a  vapor  through  a  packed  bed  as  presented  in  Sections  II 
and  III.  The  former  of  these  studies  concentrated  on  the  parameters  influencing  the  LTE 
and  two-dimensionality  of  the  transport  phenomena  while  the  latter  one  concentrates  on  the 
condensing  flows  through  a  packed  bed.  The  LTE  assumption  has  not  been  used  in  any  of 
these  studies,  and  inertia  effects  have  been  accounted  for  by  the  use  of  the  Ergur.- 
Forchheimer  relation  rather  than  Darcy  flow  formulation. 

In  the  present  section,  compressible  flow  of  an  ideal  gas  through  a  packed  bed  is 
investigated  for  oscillating  inlet  pressure  and  oscillating  inlet  temperature  boundary 
conditions  in  order  to  represent  more  closely  the  real-life  conditions  for  certain  situations. 
For  example,  in  real  applications,  more  often  than  not,  some  form  of  oscillation  prevails  in 
the  inlet  pressure  or  temperature.  It  is  then  crucial  to  know  the  qualitative  and  quantitative 
effects  of  the  oscillations  on  the  thermal  charging  characteristics  and  the  net  energy  storage 
capabilities  of  the  packed  bed.  Our  aim  is  to  analyze  the  behavior  of  the  transport  processes 
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and  energy  storage  characteristics  in  these  oscillating  boundary  condition  flows  through 
porous  media  with  specific  attention  on  packed  beds. 


6.2  Problem  statement  and  formulation 

A  rectangular  packed  bed  is  assumed  to  be  formed  by  regularly  shaped  and  sized 
spheres  packed  between  two  horizontal  walls.  The  schematic  diagram  of  the  physical 
system  under  consideration  is  depicted  in  Figure  6.1.  Relatively  high-speed  flows  are 
considered  in  the  present  study  and,  therefore,  the  flow  is  essentially  forced  convective  in 
nature.  In  the  present  study  the  top  and  bottom  walls  are  assumed  to  be  insulated,  and  the 
depth  to  be  infinitely  long,  thus,  rendering  the  problem  essentially  one-dimensional.  The 
working  fluid  was  taken  to  be  superheated  Refrigerant- 12  which  was  modeled  as  an  ideal 
gas,  while  the  material  of  the  packed  bed  particles  was  chosen  to  be  1%  carbon-steel.  The 
governing  equations  following  the  previous  section  are  given  as  follows: 

<P,>')  +  £(<fi,  >'<«,>)  =0  (6.1) 
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<Pv>v=<pv>vR<Tv>v  (6.5) 

These  equations  represent  the  gas  phase  continuity  equation,  gas  phase  momentum 
equation,  gas  phase  energy  equation,  solid  phase  energy  equation  and  the  equation  of  state 
for  the  working  fluid,  respectively,  where  the  unknown  variables  solved  from  these 
equations  are,  respectively,  <  pv  >v,  <  uv  >,  <  Tv  >v,  <  Ts  >s,  and  <  Pv  >v. 

Modeling  of  the  permeability  of  the  packed  bed,  Kv,  and  the  geometric  function  F 
appearing  in  the  gas  phase  momentum  equation,  the  effective  thermal  conductivities,  the 
fluid-to-particle  heat  transfer  coefficient  and  the  specific  interfacial  surface  area  of  the 
packed  bed  was  performed  as  in  Section  II. 


The  initial  conditions  employed  in  the  solution  of  the  problem  were: 


Pv(x,t  =  0)  =  P0 

Tv(x,t  =  0)  =  Ts(x,t  =  0)  =  T0 

uv(x,t  =  0)  =  0 


(6.6) 


corresponding  values  of  pv  being  computed  from  equation  (5). 

The  boundary  conditions  used  for  the  case  with  oscillating  inlet  gas  pressure  can  be 
mathematically  expressed  as: 


Pv(x  =  0,t)  =  Pav  +  A  Cos(2jt/t) 
Pv(x  =  L,t)  =  P0 
Tv(x  =  0,t)  =  Tvlv 


where  /is  the  frequency  and  A  is  the  amplitude  of  the  cosinusoidal  variation  of  the  inlet 
pressure,  and  again  the  corresponding  py  values  are  computed  from  equation  (5). 


The  boundary  conditions  employed  for  the  case  with  oscillating  gas  phase  inlet 
temperature  condition  can  similarly  be  expressed  as: 
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(6.8) 


Pv(x  =  0,t)  =  P„ 

Pv(x  =  L,t)  =  P0 

Tv(x  =  0,t)  =  Tvtv  +  B  Cos(27c/t) 

where  /  is  the  frequency  and  B  is  the  amplitude  of  the  variation  of  the  inlet  gas 
temperature.  pv  is  computed  as  in  the  previous  case.  The  boundary  conditions  for  the  case 
with  constant  temperature  and  constant  pressure  at  the  inlet  are  similar  to  those  in  equation 
(12)  with  A  being  equal  to  zero. 

In  the  case  studies  performed,  the  following  values  have  been  used  for  the 
variables: 

P0  =  100  kPa  P,v  =  104  kPa  T0  =  300  K  TVIV=350K. 

The  nominal  particle  Reynolds  number  for  the  average  inlet  pressure  and  temperature  was 
745. 


6.4  Results  and  discussions 

The  solution  procedure  was  exactly  similar  to  those  of  the  problems  in  the  previous 
sections.  In  order  to  explore  any  differences  in  the  transport  phenomena  and  energy 
storage  characteristics  of  the  packed  bed  with  oscillating  flow  boundary  conditions  from 
those  of  constant  boundary  conditions,  we  considered  these  cases  for  different  ranges  by 
using  different  values  for  parameters  A,/,  and  B  in  equations  (6.7)  and  (6.8)  in  different 
runs. 

In  previous  sections  we  have  observed  that  for  constant  temperature  and  pressure 
boundary  conditions  at  the  inlet,  the  solution  of  the  problem  had  two  distinct  stages. 
Namely,  the  early  stage  during  which  the  pressure  evolution  in  the  packed  bed  takes  place 
very  rapidly  causing  drastic  variation  in  the  other  field  variables,  and  the  later  stage  in 
which  the  changes  in  the  field  variables  are  mostly  temperature  dependent.  The 
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distributions  of  different  field  variables  along  the  packed  bed  during  the  early  stage  are 
depicted  in  Figure  6.2  for  the  case  with  constant  inlet  temperature  and  pressure.  In  this 
Figure  the  gas  phase  density  was  nondimensionalized  with  respect  to  an  average  density 
value  calculated  from  the  equation  of  state  using  average  values  of  gas  phase  pressure  and 
temperature,  i.e.,  average  of  the  initial  and  inlet  conditions.  Likewise,  the  nominal  value  of 
the  velocity  used  for  nondimensionalizing  the  gas  phase  velocity  was  computed  from 
equation  (6.2)  based  on  the  average  density  and  the  average  pressure  gradient  along  the 
packed  bed.  These  characteristics  have  been  discussed  in  detail  in  Sections  II  and  III. 
Similar  characteristics  were  observed  in  the  present  work  with  the  exception  that  in  the  case 
of  oscillating  inlet  pressure  condition  the  pressure  distribution  within  the  packed  bed  was 
not  very  close  linearly  as  it  was  in  the  other  cases.  Rather,  it  showed  an  oscillating 
behavior  along  the  packed  bed. 

First,  temperature  variation  of  the  solid  and  fluid  phases  within  the  packed  bed  with 
respect  to  time  was  considered.  Figure  6.3  depicts  this  variation  for  the  case  with  constant 
temperature  and  constant  pressure  inlet  boundary  conditions.  For  convenience  this  case 
will  be  called  Case  I.  Figure  6.4  shows  the  temperature  variations  for  constant 
temperature,  and  oscillating  pressure  at  the  inlet  of  the  packed  bed.  This  case  will  be  called 
Case  II.  For  Case  II,  the  chosen  parameters  were  A  =  2  kPa  and  /  =  0.05Hz.  For  Case 
III,  which  is  for  constant  pressure  but  oscillating  temperature  inlet  boundary  condition,  the 
temperature  variations  are  shown  in  Figure  6.5.  For  this  case  we  had  B  =  25K  and 
/  =  0.05 Hz.  Qualitatively  Figures  6.3  and  6.4  depict  similar  behavior  for  the  solid  and 
fluid  phase  temperature  variations  with  larger  temperature  difference  between  the  two 
phases  at  the  beginning,  and  narrowing  difference  as  the  thermal  front  moves  within  the 
packed  bed.  This  behavior  is  also  valid  for  Case  III  in  the  downstream  section  of  the 
packed  bed  but  not  at  the  entrance  region  because  of  the  oscillating  temperature  inlet 
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Figure  6.2:  Variation  of  different  field  variables  during  early  stage . 


360 


2  kPa,/ 


condition.  At  the  entrance  region  the  temperature  difference  between  the  solid  and  fluid 
phases  can  increase  or  decrease  due  to  the  oscillating  inlet  fluid  temperature  and  due  to  the 
fact  that  there  will  be  a  lag  time  for  the  response  of  the  solid  phase  temperature  to  this 
oscillation. 

Figure  6.6  depicts  the  variation  of  the  gas  phase  density,  pressure  and  velocity  at 
four  different  time  levels  for  three  different  cases  for  Case  II  during  one  complete  cycle  of 
the  pressure  oscillation.  These  three  cases  involved  different  values  for  parameters  A  and / 
as  shown  in  Figure  6.6.  Cases  in  Figures  6.6(a)  and  6.6(b)  have  the  same  amplitude  but 
different  frequency  of  inlet  pressure  oscillation,  where  cases  in  Figures  6.6(b)  and  6.6(c) 
have  the  same  frequency  but  different  amplitudes  of  inlet  pressure  oscillation.  In  each  case 
the  four  time  levels  depicted  were  chosen  so  that  they  would  span  over  one  complete  cycle 
of  oscillation  (the  very  first  cycle  in  the  process).  Although  the  pressure  distribution  is 
nearly  linear  along  the  packed  bed  duiing  the  first  three  time  levels,  a  careful  examination  of 
the  figures  at  the  fourth  time  level  reveals  that  the  pressure  distribution  in  the  packed  bed 
picks  up  from  the  oscillating  inlet  condition  and  shows  an  oscillating  behavior  along  the 
packed  bed  too,  i.e.,  spatial  oscillation  in  addition  to  temporal  oscillation.  Also  as  can  be 
expected,  the  range  of  variation  of  the  gas  phase  velocity  is  higher  in  the  case  with  larger 
amplitude  in  the  inlet  pressure  oscillation.  The  comparison  of  the  range  of  velocity 
variation  in  Figure  6.6(c)  compared  to  those  in  Figures  6.6(a)  and  6.6(b)  reveals  this 
clearly.  The  variation  in  the  gas  phase  density  can  be  explained  by  the  use  of  the  equation 
of  state  and  the  variation  of  the  pressure  and  temperature  along  the  packed  bed.  At  the 
entrance  region  the  sharp  decrease  in  the  temperature  at  the  beginning  of  the  charging 
process  requires  an  increase  in  the  gas  phase  density  since  the  rate  of  decrease  in  gas  phase 
pressure  is  less  pronounced  than  that  in  temperature. 


141 


2100 


i 

j 


o 

in 


o 

<0 

CO 


o 

t- 

CM 


O 

00 


o 

OJ 


o 


0 

0 

W 

'w' 

0 

s 

•#H 

-+-> 


o 

o 

o 

o 

o 

o 

o 

o 

© 

o 

o 

o 

00 

CM 

CO 

CO 

CO 

(pj)  paao^s  /§i9U9 


143 


Figure  6.7:  Time  history  of  energy  storage  in  the  packed  bed. 


Figure  6.7  depicts  the  time  history  of  the  net  energy  storage  per  unit  width  of  the 
packed  bed  for  Case  I  and  two  cases  of  Case  II.  As  expected,  the  asymptotic  value  of  the 
total  net  energy  stored  in  each  case  is  the  same.  As  the  amplitude  of  the  oscillations  in  the 
inlet  pressure  increases,  the  oscillations  in  the  net  energy  stored  become  more  pronounced 
since  the  oscillations  in  effect  result  in  oscillations  in  the  mass  flow  rate  of  the  gas  phase 
through  the  packed  bed.  Larger  frequencies,  on  the  other  hand,  tend  to  smooth  the 
variation  in  the  net  energy  storage. 

The  variations  in  the  rate  of  heat  flow  into  and  out  of  the  packed  bed  for  Case  I  and 
Case  II  with  A  =  2kPa  and  /  =  0.05Hz  are  depicted  in  Figure  6.8.  The  oscillating 
behavior  in  this  figure  is  that  of  Case  II  and  the  smooth  variation  is  that  of  Case  I. 
Although  the  heat  flow  rates  into  and  out  of  the  packed  bed  oscillate  in  Case  II  due  to  the 
oscillation  in  the  inlet  pressure  and,  hence,  the  mass  flow  rate,  the  average  variation  of  each 
of  these  quantities  is  qualitatively  very  similar  to  those  of  Case  I,  i.e.,  the  variation  of  the 
difference  between  the  heat  flowing  into  and  out  of  the  packed  bed  has  the  same  trend  in 
both  cases. 

The  comparison  of  the  rate  of  heat  flow  into  and  out  of  the  packed  bed  for  Cases  I 
and  III  is  shown  in  Figure  6.9.  For  Case  III  in  this  figure,  B  =  25K  and  /  =0.05Hz, 
Again,  due  to  the  oscillating  inlet  temperature,  the  density  and  velocity  of  the  gas  and, 
hence,  the  mass  flow  rate  into  the  packed  bed  oscillates.  However,  this  oscillation  in 
temperature  has  a  less  pronounced  effect  on  the  heat  flow  rates  into  and  out  of  the  packed 
bed  than  the  oscillation  in  the  inlet  pressure  has  as  shown  in  Figure  6.8.  The  difference 
between  the  scales  of  Figures  6.8  and  6.9  should  be  noted.  Part  of  the  reason  for  this 
behavior  may  be  attributed  to  the  fact  that,  in  calculating  the  gas  phase  velocity,  the 
pressure  gradient  term  in  the  gas  phase  momentum  equation  is  much  more  dominant  than 
the  inertia  term  which  involves  the  gas  phase  density,  which  in  turn  varies  proportional  to 
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Figure  6.8:  Time  history  of  heat  flow  rates  into  and  out  of  the  packed  bed. 


the  temperature.  Thus,  the  former  affects  the  variation  in  the  gas  phase  velocity  more, 
causing  a  more  pronounced  oscillation  in  the  mass  flow  rate  of  the  gas  phase. 

In  Case  III  we  cannot  talk  about  the  complete  thermal  charging  of  the  packed  bed 
due  to  the  oscillating  inlet  temperature  boundary  condition.  Because  at  the  entrance  region 
of  the  packed  bed,  due  to  this  oscillation,  the  solid  temperature  also  oscillates  (not 
necessarily  at  the  same  frequency  and  amplitude,  because  the  heat  capacity  of  the  solid 
phase  is  much  larger  than  that  of  the  working  fluid  and,  therefore,  the  temperature  of  the 
solid  phase  cannot  follow  the  temperature  of  the  incoming  gas  at  the  same  frequency  and 
amplitude).  Yet,  one  can  speak  of  a  pseudo-charging  of  the  packed  bed.  Figure  6.10 
depicts  this  pseudo-charging  behavior  (representing  the  net  energy  stored  within  the  packed 
bed).  This  behavior  is  different  from  that  of  Case  II.  Because  in  Case  II,  although  the  gas 
phase  inlet  pressure  is  oscillating,  since  the  gas  phase  inlet  temperature  is  kept  constant  at 
the  highest  value  that  can  be  attained  by  the  packed  bed  particles,  there  is  continuous 
thermal  energy  storage  within  the  packed  bed  until  the  packed  bed  is  thermally  fully 
charged.  However,  in  Case  III,  due  to  the  variation  of  the  solid  temperature  in  the  entrance 
region,  there  are  alternating  energy  storage  and  removal  from  the  packed  bed. 

In  order  to  find  out  whether  the  qualitative  behavior  of  the  energy  storage 
characteristics  changes  with  particle  Reynolds  number  (Rep),  cases  with  different  nominal 
Rep  were  investigated.  Higher  particle  Reynolds  numbers  were  obtained  by  increasing  the 
mean  inlet  pressure.  The  results  are  depicted  in  Figure  6.11  for  three  different  cases  in 
which  the  amplitude  of  the  inlet  pressure  oscillations  was  the  same.  As  can  be  seen,  the 
qualitative  behavior  is  similar  in  each  case  although  higher  nominal  Reynolds  numbers, 
meaning  higher  mass  flow  rates,  cause  faster  charging  of  the  packed  bed.  The  effect  of  the 
pressure  oscillations  is  seen  most  clearly  in  the  case  with  lowest  Rep  since  the  amplitude  of 
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the  oscillations  is  largest  relative  to  the  global  pressure.difference  applied  across  the  packed 
bed  for  that  case. 

6.5  Conclusions 

The  dynamic  response  of  sensible  heat  storage  packed  beds  with  oscillating  inlet 
boundary  conditions  has  been  studied  numerically  for  forced  convective  flow  of  a 
compressible  fluid.  A  finite  difference  scheme  with  uniform  grid  size  was  employed.  It 
was  found  that  the  average  energy  storage  behavior  did  not  have  major  differences  in  the 
cases  of  constant  or  oscillating  inlet  boundary  conditions  although  the  field  variables 
showed  oscillating  behavior  in  cases  of  oscillating  boundary  conditions.  As  expected,  the 
variation  of  the  energy  storage  was  found  to  become  smoother  as  the  amplitude  of  the 
oscillation  of  the  inlet  condition  decreased  and/or  as  the  frequency  increased.  It  was  also 
observed  that  due  to  the  nature  of  the  governing  equations  the  response  of  the  field 
variables  was  more  sensitive  to  the  amplitude  of  the  fluid  pressure  variation  than  to  that  of 
the  fluid  temperature  variation  at  the  inlet  of  the  packed  bed. 
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